diff --git a/tests/SetEpsTestCasesFromScratch/f90/Makefile b/tests/SetEpsTestCasesFromScratch/f90/Makefile new file mode 100644 index 0000000000000000000000000000000000000000..e380515f8342aff74dbf1217470d8b9ac73f5b65 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/Makefile @@ -0,0 +1,71 @@ +.SUFFIXES: .f90 .F90 .o + +# fortran compiler +FC = gfortran + +# fortran compiler options +FFLAGS = -g -gdwarf-2 -fbounds-check -fcheck=all -ffpe-trap=invalid -O0 -Wall +# gfortran version 4.8 does not know -fdiagnostics-color +# safeguard for major version >= 5 +GFORTAN_VERSION_GE_5 := $(shell echo `gcc -dumpversion | cut -f1 -d. ` \>= 5 | bc) +ifeq "$(GFORTAN_VERSION_GE_5)" "1" + FFLAGS += -fdiagnostics-color=auto +endif + +# implicit rules for compiling fortran files +%.o: %.f90 ; $(FC) $(FFLAGS) -c $< +%_mod.mod: %.f90 ; $(FC) $(FFLAGS) -c $< + +all: eels + +eelsmods = quanc8_mod.mod queels_mod.mod seteps_mod.mod param_mod.mod +doeels.o: doeels.f90 $(eelsmods) + +queels_mod.mod: quanc8_mod.mod fint1_mod.mod fint2_mod.mod fint3_mod.mod +quanc8_mod.mod: fun_mod.mod +fint1_mod.mod: surlos_mod.mod param_mod.mod +fint2_mod.mod: surlos_mod.mod param_mod.mod +fint3_mod.mod: surlos_mod.mod param_mod.mod +fun_mod.mod: param_mod.mod + +getoptsubs = getopt.o text.o system.o constants.o date_and_time.o kinds.o dummy_variables.o + +getoptmods = sufr_getopt.mod sufr_text.mod sufr_system.mod sufr_constants.mod sufr_date_and_time.mod sufr_kinds.mod sufr_dummy.mod + +sufr_getopt.mod: sufr_text.mod + $(FC) $(FFLAGS) -ffree-line-length-180 -c -o getopt.o ../getopt-libs/libsufr-0.7.7/src/getopt.f90 + +sufr_text.mod: sufr_system.mod + $(FC) $(FFLAGS) -ffree-line-length-160 -c -o text.o ../getopt-libs/libsufr-0.7.7/src/text.f90 + +sufr_system.mod: sufr_constants.mod sufr_dummy.mod + $(FC) $(FFLAGS) -c -o system.o ../getopt-libs/libsufr-0.7.7/src/system.f90 + +sufr_constants.mod: sufr_date_and_time.mod + $(FC) $(FFLAGS) -ffree-line-length-256 -c -o constants.o ../getopt-libs/libsufr-0.7.7/src/constants.f90 + +sufr_date_and_time.mod: sufr_kinds.mod + $(FC) $(FFLAGS) -ffree-line-length-150 -c -o date_and_time.o ../getopt-libs/libsufr-0.7.7/src/date_and_time.f90 + +sufr_kinds.mod: + $(FC) $(FFLAGS) -c -o kinds.o ../getopt-libs/libsufr-0.7.7/src/kinds.f90 + +sufr_dummy.mod: + $(FC) $(FFLAGS) -c -o dummy_variables.o ../getopt-libs/libsufr-0.7.7/src/dummy_variables.f90 + +eelssubs = doeels.o usurlo.o quanc8.o fun.o queels.o fint1.o fint2.o fint3.o surlos.o seteps.o phint.o qrat.o param.o +eels: eels.f90 change_working_dir.o $(eelssubs) $(eelsmods) + $(FC) $(FFLAGS) -o eels eels.f90 change_working_dir.o $(eelssubs) + +get_commandline_options.o: get_commandline_options.f90 $(getoptmods) $(getoptsubs) + $(FC) $(FFLAGS) -c -o get_commandline_options.o get_commandline_options.f90 + +clean: + rm -f *.o + rm -rf *.dSYM + rm -f *.mod + rm -f *.so + rm -f eelsou + rm -f eels eels.exe + +.PHONY: all clean diff --git a/tests/SetEpsTestCasesFromScratch/f90/README b/tests/SetEpsTestCasesFromScratch/f90/README new file mode 100644 index 0000000000000000000000000000000000000000..d4839456790eb5b2cb8d560f56616f7d3432a629 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/README @@ -0,0 +1,13 @@ +Command for compilation: + +> make + +Command for execution: + +> doRun.sh + +Cleanup: + +> make clean + +KMS \ No newline at end of file diff --git a/tests/SetEpsTestCasesFromScratch/f90/calltree.txt b/tests/SetEpsTestCasesFromScratch/f90/calltree.txt new file mode 100644 index 0000000000000000000000000000000000000000..78ec920fd02b00b1906bb0509bbd3ce2ca96625d --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/calltree.txt @@ -0,0 +1,16 @@ +Calltrees of EELS + +EELS + change_working_dir + doeels + quanc8 + fun + queels + quanc8 + fun + fint1, fint2, fint3 + usurlo + surlos + phint + seteps + (extend3) diff --git a/tests/SetEpsTestCasesFromScratch/f90/change_working_dir.f90 b/tests/SetEpsTestCasesFromScratch/f90/change_working_dir.f90 new file mode 100644 index 0000000000000000000000000000000000000000..4685606b6ae636f63c5e444486c7acb2c61e7842 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/change_working_dir.f90 @@ -0,0 +1,25 @@ +subroutine change_working_dir() + +! This routine gets the first argument of the commandline and takes it +! as the path to change the working directory +! used intrinsic routines: +! iarg returns the number of commandline arguments without the program cname. +! chdir changes the directory and returns 0 on success. +! trim removes blanks from strings. + + character (len = 256) :: argument + integer :: status + + if (iargc() == 1) then + call getarg(1, argument) + status = chdir(trim(argument)) + if (status /= 0) then + write (*,*) '*** change directory failed ***' + write (*,*) 'directory tried: ', trim(argument) + write (*,*) 'error code (see: man chdir): ', status + write (*,*) 'continuing in the start directory' + end if + end if + + return +end subroutine change_working_dir diff --git a/tests/SetEpsTestCasesFromScratch/f90/doRun.sh b/tests/SetEpsTestCasesFromScratch/f90/doRun.sh new file mode 100755 index 0000000000000000000000000000000000000000..f5a8e515567c4f4d1e9bde28d4efd172018c2af5 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/doRun.sh @@ -0,0 +1,15 @@ +#!/bin/sh -v +cp ../inputFIles/eelsin001 eelsin +./eels +mv seteps.log epsLog/seteps001.log + +cp ../inputFIles/eelsin004 eelsin +./eels +mv seteps.log epsLog/seteps004.log + +cp ../inputFIles/eelsin006 eelsin +./eels +mv seteps.log epsLog/seteps006.log + +rm eelsin +rm EELSOU diff --git a/tests/SetEpsTestCasesFromScratch/f90/doeels.f90 b/tests/SetEpsTestCasesFromScratch/f90/doeels.f90 new file mode 100644 index 0000000000000000000000000000000000000000..972e962842154d1157f35bba8928a77eebb10859 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/doeels.f90 @@ -0,0 +1,261 @@ +module doeels_mod +contains +subroutine doeels (e0, theta, phia, phib, wmin, wmax, dw, comment, comment_size, & + layers, neps, nper, name, name_size, thick, epsinf, nos, osc, osc_size,& + contrl, mode, wn_array, debug, f_array, wn_array_size) + +! ****************************************************************** +! * * +! * compute the classical eels spectrum of an arbitrary plane- * +! * statified medium made from isotropic materials in specular * +! * geometry using the dielectric theory of eels. * +! * * +! * e0: impact energy (eV) * +! * theta: incidence angle (°) * +! * phia, phib: angular apertures of the elliptic detector (°) * +! * wmin, wmax, dw: energy-loss interval and step size (cm**-1) * +! * comment: comment lines * +! * comment_size: number of comment lines * +! * layers: number of layers * +! * neps: number of epsilon values * +! * nper: number of periodic layers * +! * name: layer names * +! * name_size: number of layers * +! * thick: layer thickness * +! * epsinf: epsilon at infinite frequency * +! * nos: number of oscillators * +! * osc: oscillator parameters, wTO, Q, lambda * +! * osc_size: number of oscillators * +! * contrl: 'image' for image-charge screening * +! * mode: 'kurosawa' for kurosawa model * +! * wn_array: frequencies * +! * debug: logical * +! * f_array: spectrum * +! * wn_array_size: number of frequencies * +! * * +! ****************************************************************** + + use quanc8_mod + use queels_mod + use seteps_mod + use param_mod + + implicit none + + integer, parameter :: nt = 5 + + double precision, intent(in) :: e0, theta, phia, phib, wmin, wmax, dw + character (len = 72) :: comment(comment_size) + character (len = 10) :: name(name_size) + double precision, intent(in) :: thick(name_size), epsinf(name_size), osc(3, osc_size) + character (len = 10) :: contrl, mode + integer, intent(in) :: comment_size, name_size, osc_size, wn_array_size + integer, intent(in out) :: layers, nper, nos(name_size) + double precision, intent(in out) :: wn_array(wn_array_size), f_array(wn_array_size) + logical, intent(in), optional :: debug + + integer :: i, iw, neps, nofu, nout, nw, lmax + double precision :: a, aerr, alpha, b, beta, & + c1, c2, dx, ener, facru, f, f0, & + f1, fpic, fun, pi, prefac, psia, psii, qrat, rerr, t, & + width, wpic, x, xmin, xmax, z, z1, z2 + double precision :: table(nt) + double complex, allocatable :: eps(:) + +! **** log modification start + integer :: j +! **** log modification end + + external fun, qrat + + data aerr / 0.0d0 /, rerr / 1.0d-06 /, f / 0.0d0 /, f1 / 0.0d0 / + +! **** log modification start + open(unit = 99, file = 'seteps.log') +! **** log modification end + + if (debug) then + write (*,*) 'doeels:' + write (*,*) 'comment: ', size(comment) + write (*,*) 'name: ', size(name) + write (*,*) 'thick: ', size(thick) + write (*,*) 'epsinf: ', size(epsinf) + write (*,*) 'osc: ', size(osc), size(osc, 1), size(osc, 2) + write (*,*) 'nos: ', size(nos) + write (*,*) 'wn_array: ', size(wn_array) + write (*,*) 'f_array: ', size(f_array) + endif + +! *** machine-dependent constants +! *** epsmac + 1.0 = epsmac , cosh(argmin) = 1.0 , tanh(argmax) = 1.0 + + pi = 4 * datan(1.0d0) + + dlimf = 0.0d0 + rational = .false. + +! *** read target specifications + + user = layers == 0 + if (user) then + + if (layers == 1) rational = .true. + if (contrl == 'image') then +! *** image-charge screening factor + if (layers == 1 .and. neps == 2) then + dlimf = dsqrt(epsinf(1) * epsinf(2)) + else + dlimf = epsinf(1) + endif + dlimf = (dlimf - 1.0d0) / (dlimf + 1.0d0) + endif + endif + +! *** initialize constants + + lmax = size(thick) + nw = size(wn_array) + if (debug) write (*,*) 'lmax: ', lmax + allocate(eps(lmax)) + if (debug) write (*,*) 'eps: ', size(eps) + nout = 1 + nw / 20 + ener = 8065 * e0 + psia = phia / 180 * pi + psii = theta / 180 * pi + cospsi = dcos(psii) + sinpsi = dsin(psii) + tanpsi = dtan(psii) + prefac = dsqrt(255500 / e0)/(137 * cospsi) + facru = psia / cospsi * dsqrt(0.2624664d0 * e0) + elleps = (1.0d0 - phia / phib) * (1.0d0 + phia / phib) + acoef = sinpsi**2 + elleps * cospsi**2 + bcoef = sinpsi * cospsi + if (dlimf > 0.0d0) then + rational = .false. + if (debug) write(*,*) ' = > electron attracted by an image charge = ', dlimf +! *** dlimf : half the length unit imposed by the image force + dlimf = 1.80d0 * dlimf/(e0 * cospsi**2) + endif + if (debug) write (*,*) 'rational: ', rational + if (rational) then + +! *** set up coefficients for the rational approximation to the integral + + if (debug) write(*,*) ' = > set up a rational approximation to the integral' + call quanc8(fun, 0.0d0, pi / 2, aerr, rerr, alpha, c1, nofu, c2, eps, thick, layers, nper) + alpha = (2 / pi)**2 * alpha + c1 = 2 / pi / dsqrt(1.0d0 - elleps) * sinpsi * alpha**2 + if (c1 > 0.99d0) then + if (debug) write(*,*) ' ===> cannot do it' + rational = .false. + else + c2 = 3 * alpha**2 / (1.0d0 - c1) + c1 = c1 * c2 + xmin = wmin / (2 * ener * psia) + xmax = wmax / (2 * ener * psia) + if (xmin <= 0.0d0) xmin = 0.0d0 + dx = dmax1(0.02d0, (xmax - xmin) / nt) + z1 = 0.0d0 + z2 = 0.0d0 + do i = 1, nt + x = xmin + i * dx + call queels(x, f, aerr, rerr, facru, eps, thick, layers, nper) + table(i) = f + f = f * (1.0d0 + alpha * x)**2 + if (dabs(c2 * f - c1) < c2 * rerr) cycle + z = (1.0d0 - f) / (c2 * f - c1) + if (z <= 0.0d0) cycle + z1 = z1 + x * z * (x**2 - z) + z2 = z2 + (x * z)**2 + enddo + if (z2 == 0.0d0) then + if (debug) write(*,*) ' ===> cannot do it' + rational = .false. + else + beta = z1 / z2 + z = 0.0d0 + do i = 1, nt + x = xmin + i * dx + z = z + (table(i) - qrat(x, alpha, beta, c1, c2))**2 + enddo + z = dsqrt(z) / nt + if (z > 5.0d-03) then + if (debug) write(*,*) ' ===> cannot do it' + rational = .false. + else + if (debug) write(*, 100) alpha, c1, c2, beta, z + endif ! z > 5.0d-03 + endif ! z2 == 0.0d0 + endif ! c1 > 0.99d0 + endif ! rational + +! *** loop over the energy losses + + if (debug) write(*, 110) + +! **** log modification start + write (99, '(i5, i5)') neps, layers + write (99, '(a)') mode + do i = 1, neps + write (99, '(a, g15.7, i5)') name(i), epsinf(i), nos(i) + do j = 1, nos(i) + write (99, '(3g15.7)') osc(1,j), osc(2,j), osc(3,j) + enddo + enddo + write (99, *) +! **** log modification end + + do iw = 1, nw + f0 = f1 + f1 = f + f = 0.0d0 + wn = wmin + (iw - 1) * dw +! if (debug) write (*,*) 'wn: ', wn + if (wn >= 0.0d0) then + if (wn /= 0.0d0) then + if (.not. user) call seteps(neps, nos, osc, epsinf, wn, name, eps, layers, mode) + + x = wn / (2 * ener * psia) + if (rational) then + f = qrat(x, alpha, beta, c1, c2) * dimag(-2 / (1.0d0 + eps(1))) + else + call queels(x, f, aerr, rerr, facru, eps, thick, layers, nper) + endif + f = prefac * f / wn + endif ! wn /= 0.0d0 + + wn_array(iw) = wn + f_array(iw) = f + +! *** localize a peak using a parabolic interpolation + + if ((iw >= 3) .and. (f1 - f0 > aerr) .and. (f1 - f > aerr)) then + a = (f1 - f0) + (f1 - f) + if (a > 4 * rerr * f1) then + b = 0.5d0 * (f1 - f0 + 3 * (f1 - f)) + t = b / a + wpic = wn - t * dw + fpic = f + 0.5d0 * b * t + width = dsqrt(8 * fpic / a) * dw + if (debug) write(*, 120) wpic, fpic, width + endif ! a > 4 * rerr * f1 + endif ! iw >= 3 ... + endif ! wn >= 0.0d0 + if (mod(iw, nout) == 0) then + if (debug) write(*, 130) 100 * iw / nw, wn, f + endif + enddo + +! **** log modification start + close (unit = 99) +! **** log modification end + + return +100 format(5x, 'alpha = ', f9.4, 4x, 'c1 = ', f9.4, 4x, 'c2 = ', f9.4, 4x, & + 'beta = ', f9.4/5x, 'accuracy = ', e9.2) +110 format(//' run (%) wn (cm**-1) pcl(wn) (cm) |', & + ' peak location amplitude width') +120 format(40x, f10.2, d12.4, f10.2) +130 format(2x, f5.1, 3x, f11.3, d14.5) +end subroutine doeels +end module doeels_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/eels.f90 b/tests/SetEpsTestCasesFromScratch/f90/eels.f90 new file mode 100644 index 0000000000000000000000000000000000000000..3eb1277e01072b48e1ee02fa199b009c24b473ad --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/eels.f90 @@ -0,0 +1,160 @@ +program eels + +! ****************************************************************** +! * * +! * compute the classical eels spectrum of an arbitrary plane- * +! * statified medium made from isotropic materials in specular * +! * geometry using the dielectric theory of eels. * +! * * +! ****************************************************************** + + use doeels_mod + + implicit none + + double precision :: e0, theta, phia, phib, wmin, wmax, dw + integer :: i, j, jos, k, l, layers, neps, nper, nw + logical :: user + character (len = 72) :: comment(2) + character (len = 10) :: contrl, mode + + double precision, allocatable :: thick(:), epsinf(:), osc(:, :) + integer, allocatable :: nos(:) + character (len = 10), allocatable :: name(:) + double precision, allocatable :: wn_array(:), f(:) + + integer :: old_size_1, old_size_2 + double precision, allocatable :: tmp_osc(:, :) + logical :: debug + integer :: ioStatus + +! integer, parameter :: name_length = 10 +! *** read spectrometer parameters + + call change_working_dir() + open(unit = 11, file = 'eelsin') +! impact energy (eV) + read(11, *) e0 +! incidence angle (°) + read(11, *) theta +! angular apertures of the elliptic detector (°) + read(11, *) phia + read(11, *) phib +! energy-loss interval and step size (cm**-1) + read(11, *) wmin + read(11, *) wmax + read(11, *) dw +! comment lines + read(11, '(a72)') (comment(k), k = 1, 2) + + write(*,*) 'program eels (September 2020)' + write(*,'(a, f6.2, a, f5.1, a, f5.2, a, f5.2, a)') & + ' e0 = ', e0, ' eV, theta = ', theta, '°, phia = ', & + phia, '°, phib = ', phib, '°' + write(*,'(a, g11.4, a, g11.4, a, g11.4, a)') & + ' energy losses from', wmin, ' to', wmax, ', step = ', dw, ' cm**-1' + write(*,*) comment(1) + write(*,*) comment(2) + + if (phia <= 0.0d0 .or. phib <= 0.0d0) then + stop '*** wrong input ***' + endif + if (e0 <= 0.0d0 .or. theta + phia >= 90.0d0) then + stop '*** bad input ***' + endif + +! *** read target specifications + + read(11, *) layers, nper, mode + user = layers == 0 + if (.not. user) then + neps = layers + if (nper == -1) then + neps = layers + 1 + nper = 1 + write(*,*) 'the substrate is a anisotropic uniaxial material' + endif + if (layers < 0 .or. nper < 1 .or. nper > layers) then + stop '*** invalid target specifications ***' + endif + write(6, 101) layers, nper + jos = 0 + allocate (name(neps), thick(neps), epsinf(neps), nos(neps)) + do l = 1, neps + if (l <= layers) then + read(11, 102) name(l), thick(l) + endif + read(11, *) epsinf(l), nos(l) + write(6, 103) + if (nos(l) <= 0) then + if (l <= layers) write(6, 104) l, name(l), thick(l), epsinf(l) + if (l > layers) write(6, 105) epsinf(l) + else + do j = 1, nos(l) + if (.not. allocated(osc)) then + allocate(osc(3, nos(l))) + else if (j == 1) then + ! call extend3(osc, nos(j)) + old_size_1 = size(osc, 1) + old_size_2 = size(osc, 2) + allocate(tmp_osc(old_size_1, old_size_2 + nos(l))) + tmp_osc(1:old_size_1, 1:old_size_2) = osc + deallocate(osc) + call move_alloc(tmp_osc, osc) + endif + jos = jos + 1 + read(11, *) (osc(k, jos), k = 1, 3) + if ((j == nos(l) / 2 + 1) .and. (nos(l) > 1)) then + write(6, 106) + write(6, 107) + endif + if (j == 1) then + if (l <= layers) then + write(6, 104) l, name(l), thick(l), epsinf(l), (osc(i, jos), i = 1, 3) + else + write(6, 105) epsinf(l), (osc(i, jos), i = 1, 3) + endif + else + write(6, 108) (osc(i, jos), i = 1, 3) + endif + enddo + endif + enddo + write(*,*) + read(11, 102, IOSTAT = ioStatus) contrl + if (ioStatus /= 0) then + contrl = '' + endif + endif + + close (unit = 11) + + nw = 1 + int((wmax - wmin) / dw) + allocate (wn_array(nw), f(nw)) + + debug = .false. + call doeels(e0, theta, phia, phib, wmin, wmax, dw, comment, size(comment), & + layers, neps, nper, name, size(name), thick, epsinf, nos, osc, size(osc, 2),& + contrl, mode, wn_array, debug, f, size(wn_array)) + + open (unit = 12, file = 'eelsou') + write (12, 207) e0, theta, phia, phib + write (12, '(a72)') comment(1) + do i = 1, nw + write (12, '(2e15.7)') wn_array(i), f(i) + enddo + close (unit = 12) + + stop + +101 format(i3, ' layer(s), nper = ', i2//' l', 2x, 'material', 7x, & + 'thickness', 5x, 'epsinf', 4x, 'wto , wp', 5x, 'q', 7x, 'gam/wto') +102 format(a10, d15.5) +103 format(1x, 72('-')) +104 format(1x, i3, 2x, a10, g15.3, f10.4, f12.4, f10.4, f9.4) +105 format(31x, f10.4, f12.4, f10.4, f9.4) +106 format(45x, 'wlo , wp', 5x, 'q', 7x, 'gam/wlo') +107 format(45x, 28('-')) +108 format(41x, f12.4, f10.4, f9.4) +207 format('e0 =', f6.2, ' theta =', f5.1, ' phia =', f5.2, ' phib =', f5.2) +end program eels diff --git a/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps001.log b/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps001.log new file mode 100644 index 0000000000000000000000000000000000000000..e2ea536d49c7355a21c3ec0f3c764080d79e5d3f --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps001.log @@ -0,0 +1,331 @@ + 1 1 +No-layers +MnO 4.950000 1 + 269.0000 16.00000 0.5000000E-01 + + 50.00000 21.52103 0.1595170 + 52.00000 21.56943 0.1668695 + 54.00000 21.62002 0.1743457 + 56.00000 21.67285 0.1819523 + 58.00000 21.72795 0.1896963 + 60.00000 21.78537 0.1975850 + 62.00000 21.84514 0.2056258 + 64.00000 21.90733 0.2138265 + 66.00000 21.97197 0.2221954 + 68.00000 22.03911 0.2307409 + 70.00000 22.10881 0.2394720 + 72.00000 22.18113 0.2483979 + 74.00000 22.25612 0.2575283 + 76.00000 22.33384 0.2668734 + 78.00000 22.41437 0.2764439 + 80.00000 22.49776 0.2862509 + 82.00000 22.58410 0.2963061 + 84.00000 22.67345 0.3066220 + 86.00000 22.76589 0.3172114 + 88.00000 22.86150 0.3280879 + 90.00000 22.96038 0.3392660 + 92.00000 23.06262 0.3507607 + 94.00000 23.16830 0.3625879 + 96.00000 23.27753 0.3747646 + 98.00000 23.39042 0.3873085 + 100.0000 23.50707 0.4002384 + 102.0000 23.62761 0.4135741 + 104.0000 23.75215 0.4273368 + 106.0000 23.88082 0.4415488 + 108.0000 24.01376 0.4562338 + 110.0000 24.15112 0.4714169 + 112.0000 24.29303 0.4871248 + 114.0000 24.43967 0.5033860 + 116.0000 24.59119 0.5202306 + 118.0000 24.74778 0.5376910 + 120.0000 24.90961 0.5558015 + 122.0000 25.07689 0.5745987 + 124.0000 25.24982 0.5941219 + 126.0000 25.42862 0.6144129 + 128.0000 25.61352 0.6355167 + 130.0000 25.80477 0.6574812 + 132.0000 26.00263 0.6803581 + 134.0000 26.20736 0.7042026 + 136.0000 26.41927 0.7290744 + 138.0000 26.63867 0.7550375 + 140.0000 26.86587 0.7821607 + 142.0000 27.10123 0.8105187 + 144.0000 27.34513 0.8401916 + 146.0000 27.59796 0.8712665 + 148.0000 27.86014 0.9038373 + 150.0000 28.13212 0.9380061 + 152.0000 28.41438 0.9738835 + 154.0000 28.70744 1.011590 + 156.0000 29.01185 1.051256 + 158.0000 29.32819 1.093025 + 160.0000 29.65711 1.137052 + 162.0000 29.99927 1.183508 + 164.0000 30.35541 1.232580 + 166.0000 30.72632 1.284473 + 168.0000 31.11284 1.339411 + 170.0000 31.51589 1.397642 + 172.0000 31.93646 1.459440 + 174.0000 32.37561 1.525108 + 176.0000 32.83451 1.594979 + 178.0000 33.31441 1.669426 + 180.0000 33.81670 1.748862 + 182.0000 34.34285 1.833748 + 184.0000 34.89450 1.924598 + 186.0000 35.47344 2.021991 + 188.0000 36.08161 2.126574 + 190.0000 36.72116 2.239078 + 192.0000 37.39445 2.360328 + 194.0000 38.10408 2.491258 + 196.0000 38.85293 2.632933 + 198.0000 39.64419 2.786563 + 200.0000 40.48138 2.953537 + 202.0000 41.36846 3.135447 + 204.0000 42.30982 3.334132 + 206.0000 43.31041 3.551719 + 208.0000 44.37575 3.790683 + 210.0000 45.51211 4.053914 + 212.0000 46.72653 4.344808 + 214.0000 48.02705 4.667370 + 216.0000 49.42280 5.026353 + 218.0000 50.92424 5.427430 + 220.0000 52.54337 5.877417 + 222.0000 54.29404 6.384554 + 224.0000 56.19229 6.958881 + 226.0000 58.25682 7.612719 + 228.0000 60.50946 8.361322 + 230.0000 62.97592 9.223738 + 232.0000 65.68653 10.22400 + 234.0000 68.67733 11.39273 + 236.0000 71.99130 12.76943 + 238.0000 75.67992 14.40565 + 240.0000 79.80501 16.36962 + 242.0000 84.44097 18.75300 + 244.0000 89.67700 21.68086 + 246.0000 95.61922 25.32691 + 248.0000 102.3911 29.93688 + 250.0000 110.1299 35.86527 + 252.0000 118.9710 43.63361 + 254.0000 129.0057 54.02312 + 256.0000 140.1712 68.21886 + 258.0000 151.9839 88.01489 + 260.0000 162.9094 116.0227 + 262.0000 168.9913 155.5193 + 264.0000 161.4903 208.5716 + 266.0000 125.8030 269.3929 + 268.0000 51.76137 314.2202 + 270.0000 -41.34966 311.9429 + 272.0000 -112.3598 264.4277 + 274.0000 -145.0721 203.6377 + 276.0000 -150.9339 151.6834 + 278.0000 -144.1920 113.2758 + 280.0000 -133.0854 86.08072 + 282.0000 -121.2878 66.84453 + 284.0000 -110.2058 53.02859 + 286.0000 -100.2706 42.89900 + 288.0000 -91.52475 35.31178 + 290.0000 -83.87047 29.51224 + 292.0000 -77.17099 24.99589 + 294.0000 -71.28999 21.41924 + 296.0000 -66.10536 18.54379 + 298.0000 -61.51244 16.20070 + 300.0000 -57.42337 14.26819 + 302.0000 -53.76490 12.65691 + 304.0000 -50.47616 11.30025 + 306.0000 -47.50643 10.14782 + 308.0000 -44.81339 9.160994 + 310.0000 -42.36153 8.309762 + 312.0000 -40.12092 7.570573 + 314.0000 -38.06623 6.924736 + 316.0000 -36.17589 6.357274 + 318.0000 -34.43146 5.856081 + 320.0000 -32.81710 5.411285 + 322.0000 -31.31912 5.014779 + 324.0000 -29.92565 4.659853 + 326.0000 -28.62636 4.340919 + 328.0000 -27.41222 4.053294 + 330.0000 -26.27528 3.793027 + 332.0000 -25.20856 3.556771 + 334.0000 -24.20587 3.341674 + 336.0000 -23.26171 3.145291 + 338.0000 -22.37121 2.965522 + 340.0000 -21.53001 2.800552 + 342.0000 -20.73419 2.648806 + 344.0000 -19.98024 2.508912 + 346.0000 -19.26501 2.379672 + 348.0000 -18.58564 2.260036 + 350.0000 -17.93955 2.149077 + 352.0000 -17.32440 2.045980 + 354.0000 -16.73806 1.950021 + 356.0000 -16.17859 1.860560 + 358.0000 -15.64422 1.777024 + 360.0000 -15.13333 1.698903 + 362.0000 -14.64444 1.625741 + 364.0000 -14.17618 1.557129 + 366.0000 -13.72729 1.492698 + 368.0000 -13.29664 1.432116 + 370.0000 -12.88314 1.375085 + 372.0000 -12.48582 1.321333 + 374.0000 -12.10376 1.270614 + 376.0000 -11.73612 1.222706 + 378.0000 -11.38212 1.177405 + 380.0000 -11.04102 1.134526 + 382.0000 -10.71216 1.093901 + 384.0000 -10.39490 1.055374 + 386.0000 -10.08865 1.018805 + 388.0000 -9.792854 0.9840638 + 390.0000 -9.507008 0.9510307 + 392.0000 -9.230627 0.9195963 + 394.0000 -8.963260 0.8896596 + 396.0000 -8.704484 0.8611274 + 398.0000 -8.453902 0.8339136 + 400.0000 -8.211142 0.8079387 + 402.0000 -7.975851 0.7831290 + 404.0000 -7.747699 0.7594162 + 406.0000 -7.526376 0.7367369 + 408.0000 -7.311588 0.7150323 + 410.0000 -7.103056 0.6942477 + 412.0000 -6.900521 0.6743320 + 414.0000 -6.703733 0.6552379 + 416.0000 -6.512459 0.6369209 + 418.0000 -6.326477 0.6193398 + 420.0000 -6.145577 0.6024559 + 422.0000 -5.969560 0.5862332 + 424.0000 -5.798236 0.5706376 + 426.0000 -5.631427 0.5556376 + 428.0000 -5.468961 0.5412035 + 430.0000 -5.310677 0.5273072 + 432.0000 -5.156420 0.5139227 + 434.0000 -5.006044 0.5010252 + 436.0000 -4.859409 0.4885916 + 438.0000 -4.716382 0.4766002 + 440.0000 -4.576835 0.4650303 + 442.0000 -4.440648 0.4538626 + 444.0000 -4.307704 0.4430788 + 446.0000 -4.177894 0.4326617 + 448.0000 -4.051111 0.4225949 + 450.0000 -3.927255 0.4128631 + 452.0000 -3.806230 0.4034517 + 454.0000 -3.687943 0.3943469 + 456.0000 -3.572305 0.3855357 + 458.0000 -3.459231 0.3770056 + 460.0000 -3.348642 0.3687451 + 462.0000 -3.240458 0.3607429 + 464.0000 -3.134605 0.3529885 + 466.0000 -3.031012 0.3454718 + 468.0000 -2.929610 0.3381834 + 470.0000 -2.830333 0.3311142 + 472.0000 -2.733117 0.3242556 + 474.0000 -2.637901 0.3175994 + 476.0000 -2.544628 0.3111377 + 478.0000 -2.453240 0.3048633 + 480.0000 -2.363683 0.2987689 + 482.0000 -2.275905 0.2928479 + 484.0000 -2.189857 0.2870937 + 486.0000 -2.105488 0.2815004 + 488.0000 -2.022754 0.2760619 + 490.0000 -1.941608 0.2707727 + 492.0000 -1.862008 0.2656275 + 494.0000 -1.783912 0.2606211 + 496.0000 -1.707279 0.2557487 + 498.0000 -1.632071 0.2510056 + 500.0000 -1.558250 0.2463872 + 502.0000 -1.485780 0.2418894 + 504.0000 -1.414625 0.2375080 + 506.0000 -1.344752 0.2332391 + 508.0000 -1.276129 0.2290789 + 510.0000 -1.208723 0.2250239 + 512.0000 -1.142504 0.2210705 + 514.0000 -1.077442 0.2172154 + 516.0000 -1.013509 0.2134555 + 518.0000 -0.9506771 0.2097876 + 520.0000 -0.8889195 0.2062089 + 522.0000 -0.8282101 0.2027165 + 524.0000 -0.7685239 0.1993077 + 526.0000 -0.7098365 0.1959800 + 528.0000 -0.6521243 0.1927307 + 530.0000 -0.5953644 0.1895575 + 532.0000 -0.5395346 0.1864581 + 534.0000 -0.4846134 0.1834302 + 536.0000 -0.4305799 0.1804717 + 538.0000 -0.3774140 0.1775805 + 540.0000 -0.3250959 0.1747546 + 542.0000 -0.2736067 0.1719921 + 544.0000 -0.2229278 0.1692912 + 546.0000 -0.1730412 0.1666500 + 548.0000 -0.1239296 0.1640669 + 550.0000 -0.7557588E-01 0.1615402 + 552.0000 -0.2796367E-01 0.1590682 + 554.0000 0.1892304E-01 0.1566495 + 556.0000 0.6509980E-01 0.1542826 + 558.0000 0.1105817 0.1519659 + 560.0000 0.1553836 0.1496982 + 562.0000 0.1995196 0.1474781 + 564.0000 0.2430036 0.1453042 + 566.0000 0.2858493 0.1431754 + 568.0000 0.3280698 0.1410904 + 570.0000 0.3696778 0.1390480 + 572.0000 0.4106859 0.1370471 + 574.0000 0.4511062 0.1350867 + 576.0000 0.4909504 0.1331656 + 578.0000 0.5302303 0.1312828 + 580.0000 0.5689568 0.1294374 + 582.0000 0.6071411 0.1276283 + 584.0000 0.6447936 0.1258547 + 586.0000 0.6819248 0.1241157 + 588.0000 0.7185448 0.1224103 + 590.0000 0.7546635 0.1207377 + 592.0000 0.7902905 0.1190971 + 594.0000 0.8254350 0.1174878 + 596.0000 0.8601064 0.1159089 + 598.0000 0.8943135 0.1143597 + 600.0000 0.9280649 0.1128394 + 602.0000 0.9613692 0.1113474 + 604.0000 0.9942347 0.1098830 + 606.0000 1.026669 0.1084456 + 608.0000 1.058681 0.1070344 + 610.0000 1.090278 0.1056489 + 612.0000 1.121467 0.1042884 + 614.0000 1.152256 0.1029524 + 616.0000 1.182651 0.1016403 + 618.0000 1.212661 0.1003516 + 620.0000 1.242291 0.9908564E-01 + 622.0000 1.271549 0.9784198E-01 + 624.0000 1.300441 0.9662010E-01 + 626.0000 1.328973 0.9541951E-01 + 628.0000 1.357151 0.9423972E-01 + 630.0000 1.384982 0.9308027E-01 + 632.0000 1.412472 0.9194070E-01 + 634.0000 1.439627 0.9082058E-01 + 636.0000 1.466451 0.8971946E-01 + 638.0000 1.492951 0.8863695E-01 + 640.0000 1.519133 0.8757262E-01 + 642.0000 1.545001 0.8652609E-01 + 644.0000 1.570561 0.8549697E-01 + 646.0000 1.595818 0.8448488E-01 + 648.0000 1.620776 0.8348947E-01 + 650.0000 1.645442 0.8251037E-01 + 652.0000 1.669819 0.8154724E-01 + 654.0000 1.693912 0.8059975E-01 + 656.0000 1.717727 0.7966756E-01 + 658.0000 1.741266 0.7875036E-01 + 660.0000 1.764535 0.7784783E-01 + 662.0000 1.787539 0.7695968E-01 + 664.0000 1.810281 0.7608560E-01 + 666.0000 1.832765 0.7522531E-01 + 668.0000 1.854995 0.7437853E-01 + 670.0000 1.876976 0.7354498E-01 + 672.0000 1.898712 0.7272440E-01 + 674.0000 1.920205 0.7191652E-01 + 676.0000 1.941461 0.7112110E-01 + 678.0000 1.962482 0.7033788E-01 + 680.0000 1.983272 0.6956663E-01 + 682.0000 2.003834 0.6880710E-01 + 684.0000 2.024173 0.6805907E-01 + 686.0000 2.044291 0.6732231E-01 + 688.0000 2.064193 0.6659661E-01 + 690.0000 2.083880 0.6588175E-01 + 692.0000 2.103356 0.6517751E-01 + 694.0000 2.122625 0.6448371E-01 + 696.0000 2.141689 0.6380014E-01 + 698.0000 2.160551 0.6312659E-01 + 700.0000 2.179215 0.6246290E-01 diff --git a/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps004.log b/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps004.log new file mode 100644 index 0000000000000000000000000000000000000000..73633ac7b3ac9d63ca9ecab6c95b9469f3c1817c --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps004.log @@ -0,0 +1,331 @@ + 1 1 +No-layers +Pt 8.900000 1 + 185541.0 -1.000000 0.1200000 + + 50.00000 -60.54409 30923.34 + 52.00000 -60.54407 29733.97 + 54.00000 -60.54404 28632.70 + 56.00000 -60.54401 27610.09 + 58.00000 -60.54397 26658.01 + 60.00000 -60.54394 25769.40 + 62.00000 -60.54391 24938.11 + 64.00000 -60.54387 24158.78 + 66.00000 -60.54383 23426.69 + 68.00000 -60.54380 22737.66 + 70.00000 -60.54376 22088.00 + 72.00000 -60.54372 21474.43 + 74.00000 -60.54368 20894.03 + 76.00000 -60.54364 20344.17 + 78.00000 -60.54359 19822.51 + 80.00000 -60.54355 19326.94 + 82.00000 -60.54350 18855.54 + 84.00000 -60.54346 18406.58 + 86.00000 -60.54341 17978.51 + 88.00000 -60.54336 17569.90 + 90.00000 -60.54331 17179.44 + 92.00000 -60.54326 16805.96 + 94.00000 -60.54321 16448.38 + 96.00000 -60.54315 16105.69 + 98.00000 -60.54310 15776.99 + 100.0000 -60.54304 15461.44 + 102.0000 -60.54299 15158.26 + 104.0000 -60.54293 14866.74 + 106.0000 -60.54287 14586.23 + 108.0000 -60.54281 14316.10 + 110.0000 -60.54275 14055.79 + 112.0000 -60.54269 13804.78 + 114.0000 -60.54262 13562.58 + 116.0000 -60.54256 13328.73 + 118.0000 -60.54249 13102.81 + 120.0000 -60.54243 12884.42 + 122.0000 -60.54236 12673.19 + 124.0000 -60.54229 12468.77 + 126.0000 -60.54222 12270.84 + 128.0000 -60.54215 12079.09 + 130.0000 -60.54208 11893.25 + 132.0000 -60.54200 11713.04 + 134.0000 -60.54193 11538.20 + 136.0000 -60.54185 11368.51 + 138.0000 -60.54178 11203.74 + 140.0000 -60.54170 11043.67 + 142.0000 -60.54162 10888.11 + 144.0000 -60.54154 10736.88 + 146.0000 -60.54146 10589.78 + 148.0000 -60.54138 10446.67 + 150.0000 -60.54129 10307.37 + 152.0000 -60.54121 10171.73 + 154.0000 -60.54112 10039.62 + 156.0000 -60.54104 9910.892 + 158.0000 -60.54095 9785.425 + 160.0000 -60.54086 9663.095 + 162.0000 -60.54077 9543.785 + 164.0000 -60.54068 9427.385 + 166.0000 -60.54058 9313.790 + 168.0000 -60.54049 9202.899 + 170.0000 -60.54040 9094.617 + 172.0000 -60.54030 8988.853 + 174.0000 -60.54020 8885.521 + 176.0000 -60.54011 8784.536 + 178.0000 -60.54001 8685.821 + 180.0000 -60.53991 8589.300 + 182.0000 -60.53980 8494.899 + 184.0000 -60.53970 8402.551 + 186.0000 -60.53960 8312.189 + 188.0000 -60.53949 8223.749 + 190.0000 -60.53939 8137.171 + 192.0000 -60.53928 8052.396 + 194.0000 -60.53917 7969.369 + 196.0000 -60.53906 7888.037 + 198.0000 -60.53895 7808.347 + 200.0000 -60.53884 7730.251 + 202.0000 -60.53873 7653.702 + 204.0000 -60.53862 7578.653 + 206.0000 -60.53850 7505.061 + 208.0000 -60.53838 7432.885 + 210.0000 -60.53827 7362.083 + 212.0000 -60.53815 7292.617 + 214.0000 -60.53803 7224.449 + 216.0000 -60.53791 7157.544 + 218.0000 -60.53779 7091.866 + 220.0000 -60.53766 7027.382 + 222.0000 -60.53754 6964.060 + 224.0000 -60.53742 6901.868 + 226.0000 -60.53729 6840.777 + 228.0000 -60.53716 6780.758 + 230.0000 -60.53703 6721.783 + 232.0000 -60.53691 6663.824 + 234.0000 -60.53677 6606.856 + 236.0000 -60.53664 6550.853 + 238.0000 -60.53651 6495.791 + 240.0000 -60.53638 6441.647 + 242.0000 -60.53624 6388.398 + 244.0000 -60.53611 6336.022 + 246.0000 -60.53597 6284.497 + 248.0000 -60.53583 6233.803 + 250.0000 -60.53569 6183.920 + 252.0000 -60.53555 6134.829 + 254.0000 -60.53541 6086.511 + 256.0000 -60.53526 6038.948 + 258.0000 -60.53512 5992.122 + 260.0000 -60.53498 5946.016 + 262.0000 -60.53483 5900.614 + 264.0000 -60.53468 5855.900 + 266.0000 -60.53453 5811.858 + 268.0000 -60.53438 5768.474 + 270.0000 -60.53423 5725.732 + 272.0000 -60.53408 5683.619 + 274.0000 -60.53393 5642.120 + 276.0000 -60.53377 5601.223 + 278.0000 -60.53362 5560.914 + 280.0000 -60.53346 5521.180 + 282.0000 -60.53331 5482.011 + 284.0000 -60.53315 5443.393 + 286.0000 -60.53299 5405.314 + 288.0000 -60.53283 5367.765 + 290.0000 -60.53267 5330.734 + 292.0000 -60.53250 5294.209 + 294.0000 -60.53234 5258.182 + 296.0000 -60.53217 5222.641 + 298.0000 -60.53201 5187.577 + 300.0000 -60.53184 5152.981 + 302.0000 -60.53167 5118.843 + 304.0000 -60.53150 5085.154 + 306.0000 -60.53133 5051.905 + 308.0000 -60.53116 5019.088 + 310.0000 -60.53098 4986.695 + 312.0000 -60.53081 4954.716 + 314.0000 -60.53064 4923.145 + 316.0000 -60.53046 4891.973 + 318.0000 -60.53028 4861.194 + 320.0000 -60.53010 4830.799 + 322.0000 -60.52992 4800.782 + 324.0000 -60.52974 4771.135 + 326.0000 -60.52956 4741.852 + 328.0000 -60.52938 4712.925 + 330.0000 -60.52919 4684.350 + 332.0000 -60.52901 4656.118 + 334.0000 -60.52882 4628.225 + 336.0000 -60.52863 4600.664 + 338.0000 -60.52844 4573.428 + 340.0000 -60.52825 4546.513 + 342.0000 -60.52806 4519.913 + 344.0000 -60.52787 4493.622 + 346.0000 -60.52768 4467.635 + 348.0000 -60.52748 4441.946 + 350.0000 -60.52729 4416.551 + 352.0000 -60.52709 4391.445 + 354.0000 -60.52689 4366.622 + 356.0000 -60.52670 4342.078 + 358.0000 -60.52650 4317.808 + 360.0000 -60.52629 4293.808 + 362.0000 -60.52609 4270.073 + 364.0000 -60.52589 4246.599 + 366.0000 -60.52568 4223.381 + 368.0000 -60.52548 4200.415 + 370.0000 -60.52527 4177.698 + 372.0000 -60.52506 4155.224 + 374.0000 -60.52486 4132.992 + 376.0000 -60.52465 4110.995 + 378.0000 -60.52443 4089.231 + 380.0000 -60.52422 4067.697 + 382.0000 -60.52401 4046.387 + 384.0000 -60.52379 4025.300 + 386.0000 -60.52358 4004.431 + 388.0000 -60.52336 3983.777 + 390.0000 -60.52314 3963.335 + 392.0000 -60.52292 3943.102 + 394.0000 -60.52270 3923.074 + 396.0000 -60.52248 3903.248 + 398.0000 -60.52226 3883.621 + 400.0000 -60.52204 3864.190 + 402.0000 -60.52181 3844.953 + 404.0000 -60.52159 3825.906 + 406.0000 -60.52136 3807.047 + 408.0000 -60.52113 3788.372 + 410.0000 -60.52090 3769.880 + 412.0000 -60.52067 3751.567 + 414.0000 -60.52044 3733.431 + 416.0000 -60.52021 3715.470 + 418.0000 -60.51998 3697.680 + 420.0000 -60.51974 3680.060 + 422.0000 -60.51951 3662.606 + 424.0000 -60.51927 3645.317 + 426.0000 -60.51903 3628.191 + 428.0000 -60.51879 3611.224 + 430.0000 -60.51855 3594.415 + 432.0000 -60.51831 3577.762 + 434.0000 -60.51807 3561.262 + 436.0000 -60.51782 3544.914 + 438.0000 -60.51758 3528.714 + 440.0000 -60.51733 3512.662 + 442.0000 -60.51709 3496.755 + 444.0000 -60.51684 3480.992 + 446.0000 -60.51659 3465.370 + 448.0000 -60.51634 3449.887 + 450.0000 -60.51609 3434.541 + 452.0000 -60.51584 3419.332 + 454.0000 -60.51558 3404.256 + 456.0000 -60.51533 3389.313 + 458.0000 -60.51507 3374.500 + 460.0000 -60.51481 3359.816 + 462.0000 -60.51456 3345.259 + 464.0000 -60.51430 3330.827 + 466.0000 -60.51404 3316.519 + 468.0000 -60.51378 3302.334 + 470.0000 -60.51351 3288.269 + 472.0000 -60.51325 3274.323 + 474.0000 -60.51298 3260.495 + 476.0000 -60.51272 3246.783 + 478.0000 -60.51245 3233.186 + 480.0000 -60.51218 3219.701 + 482.0000 -60.51191 3206.329 + 484.0000 -60.51164 3193.068 + 486.0000 -60.51137 3179.915 + 488.0000 -60.51110 3166.870 + 490.0000 -60.51083 3153.932 + 492.0000 -60.51055 3141.098 + 494.0000 -60.51028 3128.369 + 496.0000 -60.51000 3115.742 + 498.0000 -60.50972 3103.217 + 500.0000 -60.50944 3090.791 + 502.0000 -60.50916 3078.465 + 504.0000 -60.50888 3066.236 + 506.0000 -60.50860 3054.104 + 508.0000 -60.50831 3042.068 + 510.0000 -60.50803 3030.126 + 512.0000 -60.50774 3018.277 + 514.0000 -60.50745 3006.520 + 516.0000 -60.50717 2994.855 + 518.0000 -60.50688 2983.279 + 520.0000 -60.50659 2971.792 + 522.0000 -60.50629 2960.394 + 524.0000 -60.50600 2949.082 + 526.0000 -60.50571 2937.857 + 528.0000 -60.50541 2926.716 + 530.0000 -60.50512 2915.659 + 532.0000 -60.50482 2904.686 + 534.0000 -60.50452 2893.794 + 536.0000 -60.50422 2882.984 + 538.0000 -60.50392 2872.254 + 540.0000 -60.50362 2861.604 + 542.0000 -60.50332 2851.032 + 544.0000 -60.50301 2840.538 + 546.0000 -60.50271 2830.120 + 548.0000 -60.50240 2819.779 + 550.0000 -60.50209 2809.513 + 552.0000 -60.50179 2799.321 + 554.0000 -60.50148 2789.203 + 556.0000 -60.50117 2779.157 + 558.0000 -60.50085 2769.184 + 560.0000 -60.50054 2759.281 + 562.0000 -60.50023 2749.449 + 564.0000 -60.49991 2739.687 + 566.0000 -60.49960 2729.994 + 568.0000 -60.49928 2720.369 + 570.0000 -60.49896 2710.811 + 572.0000 -60.49864 2701.320 + 574.0000 -60.49832 2691.896 + 576.0000 -60.49800 2682.536 + 578.0000 -60.49768 2673.242 + 580.0000 -60.49735 2664.011 + 582.0000 -60.49703 2654.844 + 584.0000 -60.49670 2645.740 + 586.0000 -60.49637 2636.697 + 588.0000 -60.49604 2627.717 + 590.0000 -60.49571 2618.797 + 592.0000 -60.49538 2609.937 + 594.0000 -60.49505 2601.137 + 596.0000 -60.49472 2592.396 + 598.0000 -60.49439 2583.713 + 600.0000 -60.49405 2575.088 + 602.0000 -60.49371 2566.521 + 604.0000 -60.49338 2558.010 + 606.0000 -60.49304 2549.555 + 608.0000 -60.49270 2541.156 + 610.0000 -60.49236 2532.812 + 612.0000 -60.49202 2524.522 + 614.0000 -60.49167 2516.287 + 616.0000 -60.49133 2508.105 + 618.0000 -60.49098 2499.975 + 620.0000 -60.49064 2491.898 + 622.0000 -60.49029 2483.873 + 624.0000 -60.48994 2475.900 + 626.0000 -60.48959 2467.977 + 628.0000 -60.48924 2460.105 + 630.0000 -60.48889 2452.283 + 632.0000 -60.48854 2444.510 + 634.0000 -60.48818 2436.786 + 636.0000 -60.48783 2429.111 + 638.0000 -60.48747 2421.483 + 640.0000 -60.48711 2413.904 + 642.0000 -60.48675 2406.372 + 644.0000 -60.48639 2398.886 + 646.0000 -60.48603 2391.447 + 648.0000 -60.48567 2384.053 + 650.0000 -60.48531 2376.705 + 652.0000 -60.48494 2369.402 + 654.0000 -60.48458 2362.144 + 656.0000 -60.48421 2354.930 + 658.0000 -60.48385 2347.760 + 660.0000 -60.48348 2340.633 + 662.0000 -60.48311 2333.549 + 664.0000 -60.48274 2326.508 + 666.0000 -60.48236 2319.509 + 668.0000 -60.48199 2312.552 + 670.0000 -60.48162 2305.636 + 672.0000 -60.48124 2298.762 + 674.0000 -60.48086 2291.928 + 676.0000 -60.48049 2285.135 + 678.0000 -60.48011 2278.381 + 680.0000 -60.47973 2271.668 + 682.0000 -60.47935 2264.994 + 684.0000 -60.47897 2258.358 + 686.0000 -60.47858 2251.762 + 688.0000 -60.47820 2245.204 + 690.0000 -60.47781 2238.683 + 692.0000 -60.47743 2232.201 + 694.0000 -60.47704 2225.755 + 696.0000 -60.47665 2219.347 + 698.0000 -60.47626 2212.976 + 700.0000 -60.47587 2206.640 diff --git a/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps006.log b/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps006.log new file mode 100644 index 0000000000000000000000000000000000000000..f050380713a3a3608f92367fb8090e3595a5f56c --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/epsLog/seteps006.log @@ -0,0 +1,333 @@ + 2 2 +No-layers +MnO 4.950000 1 + 269.0000 16.00000 0.5000000E-01 +Platinum 8.900000 1 + 269.0000 16.00000 0.5000000E-01 + + 50.00000 21.52103 0.1595170 -60.54409 30923.34 + 52.00000 21.56943 0.1668695 -60.54407 29733.97 + 54.00000 21.62002 0.1743457 -60.54404 28632.70 + 56.00000 21.67285 0.1819523 -60.54401 27610.09 + 58.00000 21.72795 0.1896963 -60.54397 26658.01 + 60.00000 21.78537 0.1975850 -60.54394 25769.40 + 62.00000 21.84514 0.2056258 -60.54391 24938.11 + 64.00000 21.90733 0.2138265 -60.54387 24158.78 + 66.00000 21.97197 0.2221954 -60.54383 23426.69 + 68.00000 22.03911 0.2307409 -60.54380 22737.66 + 70.00000 22.10881 0.2394720 -60.54376 22088.00 + 72.00000 22.18113 0.2483979 -60.54372 21474.43 + 74.00000 22.25612 0.2575283 -60.54368 20894.03 + 76.00000 22.33384 0.2668734 -60.54364 20344.17 + 78.00000 22.41437 0.2764439 -60.54359 19822.51 + 80.00000 22.49776 0.2862509 -60.54355 19326.94 + 82.00000 22.58410 0.2963061 -60.54350 18855.54 + 84.00000 22.67345 0.3066220 -60.54346 18406.58 + 86.00000 22.76589 0.3172114 -60.54341 17978.51 + 88.00000 22.86150 0.3280879 -60.54336 17569.90 + 90.00000 22.96038 0.3392660 -60.54331 17179.44 + 92.00000 23.06262 0.3507607 -60.54326 16805.96 + 94.00000 23.16830 0.3625879 -60.54321 16448.38 + 96.00000 23.27753 0.3747646 -60.54315 16105.69 + 98.00000 23.39042 0.3873085 -60.54310 15776.99 + 100.0000 23.50707 0.4002384 -60.54304 15461.44 + 102.0000 23.62761 0.4135741 -60.54299 15158.26 + 104.0000 23.75215 0.4273368 -60.54293 14866.74 + 106.0000 23.88082 0.4415488 -60.54287 14586.23 + 108.0000 24.01376 0.4562338 -60.54281 14316.10 + 110.0000 24.15112 0.4714169 -60.54275 14055.79 + 112.0000 24.29303 0.4871248 -60.54269 13804.78 + 114.0000 24.43967 0.5033860 -60.54262 13562.58 + 116.0000 24.59119 0.5202306 -60.54256 13328.73 + 118.0000 24.74778 0.5376910 -60.54249 13102.81 + 120.0000 24.90961 0.5558015 -60.54243 12884.42 + 122.0000 25.07689 0.5745987 -60.54236 12673.19 + 124.0000 25.24982 0.5941219 -60.54229 12468.77 + 126.0000 25.42862 0.6144129 -60.54222 12270.84 + 128.0000 25.61352 0.6355167 -60.54215 12079.09 + 130.0000 25.80477 0.6574812 -60.54208 11893.25 + 132.0000 26.00263 0.6803581 -60.54200 11713.04 + 134.0000 26.20736 0.7042026 -60.54193 11538.20 + 136.0000 26.41927 0.7290744 -60.54185 11368.51 + 138.0000 26.63867 0.7550375 -60.54178 11203.74 + 140.0000 26.86587 0.7821607 -60.54170 11043.67 + 142.0000 27.10123 0.8105187 -60.54162 10888.11 + 144.0000 27.34513 0.8401916 -60.54154 10736.88 + 146.0000 27.59796 0.8712665 -60.54146 10589.78 + 148.0000 27.86014 0.9038373 -60.54138 10446.67 + 150.0000 28.13212 0.9380061 -60.54129 10307.37 + 152.0000 28.41438 0.9738835 -60.54121 10171.73 + 154.0000 28.70744 1.011590 -60.54112 10039.62 + 156.0000 29.01185 1.051256 -60.54104 9910.892 + 158.0000 29.32819 1.093025 -60.54095 9785.425 + 160.0000 29.65711 1.137052 -60.54086 9663.095 + 162.0000 29.99927 1.183508 -60.54077 9543.785 + 164.0000 30.35541 1.232580 -60.54068 9427.385 + 166.0000 30.72632 1.284473 -60.54058 9313.790 + 168.0000 31.11284 1.339411 -60.54049 9202.899 + 170.0000 31.51589 1.397642 -60.54040 9094.617 + 172.0000 31.93646 1.459440 -60.54030 8988.853 + 174.0000 32.37561 1.525108 -60.54020 8885.521 + 176.0000 32.83451 1.594979 -60.54011 8784.536 + 178.0000 33.31441 1.669426 -60.54001 8685.821 + 180.0000 33.81670 1.748862 -60.53991 8589.300 + 182.0000 34.34285 1.833748 -60.53980 8494.899 + 184.0000 34.89450 1.924598 -60.53970 8402.551 + 186.0000 35.47344 2.021991 -60.53960 8312.189 + 188.0000 36.08161 2.126574 -60.53949 8223.749 + 190.0000 36.72116 2.239078 -60.53939 8137.171 + 192.0000 37.39445 2.360328 -60.53928 8052.396 + 194.0000 38.10408 2.491258 -60.53917 7969.369 + 196.0000 38.85293 2.632933 -60.53906 7888.037 + 198.0000 39.64419 2.786563 -60.53895 7808.347 + 200.0000 40.48138 2.953537 -60.53884 7730.251 + 202.0000 41.36846 3.135447 -60.53873 7653.702 + 204.0000 42.30982 3.334132 -60.53862 7578.653 + 206.0000 43.31041 3.551719 -60.53850 7505.061 + 208.0000 44.37575 3.790683 -60.53838 7432.885 + 210.0000 45.51211 4.053914 -60.53827 7362.083 + 212.0000 46.72653 4.344808 -60.53815 7292.617 + 214.0000 48.02705 4.667370 -60.53803 7224.449 + 216.0000 49.42280 5.026353 -60.53791 7157.544 + 218.0000 50.92424 5.427430 -60.53779 7091.866 + 220.0000 52.54337 5.877417 -60.53766 7027.382 + 222.0000 54.29404 6.384554 -60.53754 6964.060 + 224.0000 56.19229 6.958881 -60.53742 6901.868 + 226.0000 58.25682 7.612719 -60.53729 6840.777 + 228.0000 60.50946 8.361322 -60.53716 6780.758 + 230.0000 62.97592 9.223738 -60.53703 6721.783 + 232.0000 65.68653 10.22400 -60.53691 6663.824 + 234.0000 68.67733 11.39273 -60.53677 6606.856 + 236.0000 71.99130 12.76943 -60.53664 6550.853 + 238.0000 75.67992 14.40565 -60.53651 6495.791 + 240.0000 79.80501 16.36962 -60.53638 6441.647 + 242.0000 84.44097 18.75300 -60.53624 6388.398 + 244.0000 89.67700 21.68086 -60.53611 6336.022 + 246.0000 95.61922 25.32691 -60.53597 6284.497 + 248.0000 102.3911 29.93688 -60.53583 6233.803 + 250.0000 110.1299 35.86527 -60.53569 6183.920 + 252.0000 118.9710 43.63361 -60.53555 6134.829 + 254.0000 129.0057 54.02312 -60.53541 6086.511 + 256.0000 140.1712 68.21886 -60.53526 6038.948 + 258.0000 151.9839 88.01489 -60.53512 5992.122 + 260.0000 162.9094 116.0227 -60.53498 5946.016 + 262.0000 168.9913 155.5193 -60.53483 5900.614 + 264.0000 161.4903 208.5716 -60.53468 5855.900 + 266.0000 125.8030 269.3929 -60.53453 5811.858 + 268.0000 51.76137 314.2202 -60.53438 5768.474 + 270.0000 -41.34966 311.9429 -60.53423 5725.732 + 272.0000 -112.3598 264.4277 -60.53408 5683.619 + 274.0000 -145.0721 203.6377 -60.53393 5642.120 + 276.0000 -150.9339 151.6834 -60.53377 5601.223 + 278.0000 -144.1920 113.2758 -60.53362 5560.914 + 280.0000 -133.0854 86.08072 -60.53346 5521.180 + 282.0000 -121.2878 66.84453 -60.53331 5482.011 + 284.0000 -110.2058 53.02859 -60.53315 5443.393 + 286.0000 -100.2706 42.89900 -60.53299 5405.314 + 288.0000 -91.52475 35.31178 -60.53283 5367.765 + 290.0000 -83.87047 29.51224 -60.53267 5330.734 + 292.0000 -77.17099 24.99589 -60.53250 5294.209 + 294.0000 -71.28999 21.41924 -60.53234 5258.182 + 296.0000 -66.10536 18.54379 -60.53217 5222.641 + 298.0000 -61.51244 16.20070 -60.53201 5187.577 + 300.0000 -57.42337 14.26819 -60.53184 5152.981 + 302.0000 -53.76490 12.65691 -60.53167 5118.843 + 304.0000 -50.47616 11.30025 -60.53150 5085.154 + 306.0000 -47.50643 10.14782 -60.53133 5051.905 + 308.0000 -44.81339 9.160994 -60.53116 5019.088 + 310.0000 -42.36153 8.309762 -60.53098 4986.695 + 312.0000 -40.12092 7.570573 -60.53081 4954.716 + 314.0000 -38.06623 6.924736 -60.53064 4923.145 + 316.0000 -36.17589 6.357274 -60.53046 4891.973 + 318.0000 -34.43146 5.856081 -60.53028 4861.194 + 320.0000 -32.81710 5.411285 -60.53010 4830.799 + 322.0000 -31.31912 5.014779 -60.52992 4800.782 + 324.0000 -29.92565 4.659853 -60.52974 4771.135 + 326.0000 -28.62636 4.340919 -60.52956 4741.852 + 328.0000 -27.41222 4.053294 -60.52938 4712.925 + 330.0000 -26.27528 3.793027 -60.52919 4684.350 + 332.0000 -25.20856 3.556771 -60.52901 4656.118 + 334.0000 -24.20587 3.341674 -60.52882 4628.225 + 336.0000 -23.26171 3.145291 -60.52863 4600.664 + 338.0000 -22.37121 2.965522 -60.52844 4573.428 + 340.0000 -21.53001 2.800552 -60.52825 4546.513 + 342.0000 -20.73419 2.648806 -60.52806 4519.913 + 344.0000 -19.98024 2.508912 -60.52787 4493.622 + 346.0000 -19.26501 2.379672 -60.52768 4467.635 + 348.0000 -18.58564 2.260036 -60.52748 4441.946 + 350.0000 -17.93955 2.149077 -60.52729 4416.551 + 352.0000 -17.32440 2.045980 -60.52709 4391.445 + 354.0000 -16.73806 1.950021 -60.52689 4366.622 + 356.0000 -16.17859 1.860560 -60.52670 4342.078 + 358.0000 -15.64422 1.777024 -60.52650 4317.808 + 360.0000 -15.13333 1.698903 -60.52629 4293.808 + 362.0000 -14.64444 1.625741 -60.52609 4270.073 + 364.0000 -14.17618 1.557129 -60.52589 4246.599 + 366.0000 -13.72729 1.492698 -60.52568 4223.381 + 368.0000 -13.29664 1.432116 -60.52548 4200.415 + 370.0000 -12.88314 1.375085 -60.52527 4177.698 + 372.0000 -12.48582 1.321333 -60.52506 4155.224 + 374.0000 -12.10376 1.270614 -60.52486 4132.992 + 376.0000 -11.73612 1.222706 -60.52465 4110.995 + 378.0000 -11.38212 1.177405 -60.52443 4089.231 + 380.0000 -11.04102 1.134526 -60.52422 4067.697 + 382.0000 -10.71216 1.093901 -60.52401 4046.387 + 384.0000 -10.39490 1.055374 -60.52379 4025.300 + 386.0000 -10.08865 1.018805 -60.52358 4004.431 + 388.0000 -9.792854 0.9840638 -60.52336 3983.777 + 390.0000 -9.507008 0.9510307 -60.52314 3963.335 + 392.0000 -9.230627 0.9195963 -60.52292 3943.102 + 394.0000 -8.963260 0.8896596 -60.52270 3923.074 + 396.0000 -8.704484 0.8611274 -60.52248 3903.248 + 398.0000 -8.453902 0.8339136 -60.52226 3883.621 + 400.0000 -8.211142 0.8079387 -60.52204 3864.190 + 402.0000 -7.975851 0.7831290 -60.52181 3844.953 + 404.0000 -7.747699 0.7594162 -60.52159 3825.906 + 406.0000 -7.526376 0.7367369 -60.52136 3807.047 + 408.0000 -7.311588 0.7150323 -60.52113 3788.372 + 410.0000 -7.103056 0.6942477 -60.52090 3769.880 + 412.0000 -6.900521 0.6743320 -60.52067 3751.567 + 414.0000 -6.703733 0.6552379 -60.52044 3733.431 + 416.0000 -6.512459 0.6369209 -60.52021 3715.470 + 418.0000 -6.326477 0.6193398 -60.51998 3697.680 + 420.0000 -6.145577 0.6024559 -60.51974 3680.060 + 422.0000 -5.969560 0.5862332 -60.51951 3662.606 + 424.0000 -5.798236 0.5706376 -60.51927 3645.317 + 426.0000 -5.631427 0.5556376 -60.51903 3628.191 + 428.0000 -5.468961 0.5412035 -60.51879 3611.224 + 430.0000 -5.310677 0.5273072 -60.51855 3594.415 + 432.0000 -5.156420 0.5139227 -60.51831 3577.762 + 434.0000 -5.006044 0.5010252 -60.51807 3561.262 + 436.0000 -4.859409 0.4885916 -60.51782 3544.914 + 438.0000 -4.716382 0.4766002 -60.51758 3528.714 + 440.0000 -4.576835 0.4650303 -60.51733 3512.662 + 442.0000 -4.440648 0.4538626 -60.51709 3496.755 + 444.0000 -4.307704 0.4430788 -60.51684 3480.992 + 446.0000 -4.177894 0.4326617 -60.51659 3465.370 + 448.0000 -4.051111 0.4225949 -60.51634 3449.887 + 450.0000 -3.927255 0.4128631 -60.51609 3434.541 + 452.0000 -3.806230 0.4034517 -60.51584 3419.332 + 454.0000 -3.687943 0.3943469 -60.51558 3404.256 + 456.0000 -3.572305 0.3855357 -60.51533 3389.313 + 458.0000 -3.459231 0.3770056 -60.51507 3374.500 + 460.0000 -3.348642 0.3687451 -60.51481 3359.816 + 462.0000 -3.240458 0.3607429 -60.51456 3345.259 + 464.0000 -3.134605 0.3529885 -60.51430 3330.827 + 466.0000 -3.031012 0.3454718 -60.51404 3316.519 + 468.0000 -2.929610 0.3381834 -60.51378 3302.334 + 470.0000 -2.830333 0.3311142 -60.51351 3288.269 + 472.0000 -2.733117 0.3242556 -60.51325 3274.323 + 474.0000 -2.637901 0.3175994 -60.51298 3260.495 + 476.0000 -2.544628 0.3111377 -60.51272 3246.783 + 478.0000 -2.453240 0.3048633 -60.51245 3233.186 + 480.0000 -2.363683 0.2987689 -60.51218 3219.701 + 482.0000 -2.275905 0.2928479 -60.51191 3206.329 + 484.0000 -2.189857 0.2870937 -60.51164 3193.068 + 486.0000 -2.105488 0.2815004 -60.51137 3179.915 + 488.0000 -2.022754 0.2760619 -60.51110 3166.870 + 490.0000 -1.941608 0.2707727 -60.51083 3153.932 + 492.0000 -1.862008 0.2656275 -60.51055 3141.098 + 494.0000 -1.783912 0.2606211 -60.51028 3128.369 + 496.0000 -1.707279 0.2557487 -60.51000 3115.742 + 498.0000 -1.632071 0.2510056 -60.50972 3103.217 + 500.0000 -1.558250 0.2463872 -60.50944 3090.791 + 502.0000 -1.485780 0.2418894 -60.50916 3078.465 + 504.0000 -1.414625 0.2375080 -60.50888 3066.236 + 506.0000 -1.344752 0.2332391 -60.50860 3054.104 + 508.0000 -1.276129 0.2290789 -60.50831 3042.068 + 510.0000 -1.208723 0.2250239 -60.50803 3030.126 + 512.0000 -1.142504 0.2210705 -60.50774 3018.277 + 514.0000 -1.077442 0.2172154 -60.50745 3006.520 + 516.0000 -1.013509 0.2134555 -60.50717 2994.855 + 518.0000 -0.9506771 0.2097876 -60.50688 2983.279 + 520.0000 -0.8889195 0.2062089 -60.50659 2971.792 + 522.0000 -0.8282101 0.2027165 -60.50629 2960.394 + 524.0000 -0.7685239 0.1993077 -60.50600 2949.082 + 526.0000 -0.7098365 0.1959800 -60.50571 2937.857 + 528.0000 -0.6521243 0.1927307 -60.50541 2926.716 + 530.0000 -0.5953644 0.1895575 -60.50512 2915.659 + 532.0000 -0.5395346 0.1864581 -60.50482 2904.686 + 534.0000 -0.4846134 0.1834302 -60.50452 2893.794 + 536.0000 -0.4305799 0.1804717 -60.50422 2882.984 + 538.0000 -0.3774140 0.1775805 -60.50392 2872.254 + 540.0000 -0.3250959 0.1747546 -60.50362 2861.604 + 542.0000 -0.2736067 0.1719921 -60.50332 2851.032 + 544.0000 -0.2229278 0.1692912 -60.50301 2840.538 + 546.0000 -0.1730412 0.1666500 -60.50271 2830.120 + 548.0000 -0.1239296 0.1640669 -60.50240 2819.779 + 550.0000 -0.7557588E-01 0.1615402 -60.50209 2809.513 + 552.0000 -0.2796367E-01 0.1590682 -60.50179 2799.321 + 554.0000 0.1892304E-01 0.1566495 -60.50148 2789.203 + 556.0000 0.6509980E-01 0.1542826 -60.50117 2779.157 + 558.0000 0.1105817 0.1519659 -60.50085 2769.184 + 560.0000 0.1553836 0.1496982 -60.50054 2759.281 + 562.0000 0.1995196 0.1474781 -60.50023 2749.449 + 564.0000 0.2430036 0.1453042 -60.49991 2739.687 + 566.0000 0.2858493 0.1431754 -60.49960 2729.994 + 568.0000 0.3280698 0.1410904 -60.49928 2720.369 + 570.0000 0.3696778 0.1390480 -60.49896 2710.811 + 572.0000 0.4106859 0.1370471 -60.49864 2701.320 + 574.0000 0.4511062 0.1350867 -60.49832 2691.896 + 576.0000 0.4909504 0.1331656 -60.49800 2682.536 + 578.0000 0.5302303 0.1312828 -60.49768 2673.242 + 580.0000 0.5689568 0.1294374 -60.49735 2664.011 + 582.0000 0.6071411 0.1276283 -60.49703 2654.844 + 584.0000 0.6447936 0.1258547 -60.49670 2645.740 + 586.0000 0.6819248 0.1241157 -60.49637 2636.697 + 588.0000 0.7185448 0.1224103 -60.49604 2627.717 + 590.0000 0.7546635 0.1207377 -60.49571 2618.797 + 592.0000 0.7902905 0.1190971 -60.49538 2609.937 + 594.0000 0.8254350 0.1174878 -60.49505 2601.137 + 596.0000 0.8601064 0.1159089 -60.49472 2592.396 + 598.0000 0.8943135 0.1143597 -60.49439 2583.713 + 600.0000 0.9280649 0.1128394 -60.49405 2575.088 + 602.0000 0.9613692 0.1113474 -60.49371 2566.521 + 604.0000 0.9942347 0.1098830 -60.49338 2558.010 + 606.0000 1.026669 0.1084456 -60.49304 2549.555 + 608.0000 1.058681 0.1070344 -60.49270 2541.156 + 610.0000 1.090278 0.1056489 -60.49236 2532.812 + 612.0000 1.121467 0.1042884 -60.49202 2524.522 + 614.0000 1.152256 0.1029524 -60.49167 2516.287 + 616.0000 1.182651 0.1016403 -60.49133 2508.105 + 618.0000 1.212661 0.1003516 -60.49098 2499.975 + 620.0000 1.242291 0.9908564E-01 -60.49064 2491.898 + 622.0000 1.271549 0.9784198E-01 -60.49029 2483.873 + 624.0000 1.300441 0.9662010E-01 -60.48994 2475.900 + 626.0000 1.328973 0.9541951E-01 -60.48959 2467.977 + 628.0000 1.357151 0.9423972E-01 -60.48924 2460.105 + 630.0000 1.384982 0.9308027E-01 -60.48889 2452.283 + 632.0000 1.412472 0.9194070E-01 -60.48854 2444.510 + 634.0000 1.439627 0.9082058E-01 -60.48818 2436.786 + 636.0000 1.466451 0.8971946E-01 -60.48783 2429.111 + 638.0000 1.492951 0.8863695E-01 -60.48747 2421.483 + 640.0000 1.519133 0.8757262E-01 -60.48711 2413.904 + 642.0000 1.545001 0.8652609E-01 -60.48675 2406.372 + 644.0000 1.570561 0.8549697E-01 -60.48639 2398.886 + 646.0000 1.595818 0.8448488E-01 -60.48603 2391.447 + 648.0000 1.620776 0.8348947E-01 -60.48567 2384.053 + 650.0000 1.645442 0.8251037E-01 -60.48531 2376.705 + 652.0000 1.669819 0.8154724E-01 -60.48494 2369.402 + 654.0000 1.693912 0.8059975E-01 -60.48458 2362.144 + 656.0000 1.717727 0.7966756E-01 -60.48421 2354.930 + 658.0000 1.741266 0.7875036E-01 -60.48385 2347.760 + 660.0000 1.764535 0.7784783E-01 -60.48348 2340.633 + 662.0000 1.787539 0.7695968E-01 -60.48311 2333.549 + 664.0000 1.810281 0.7608560E-01 -60.48274 2326.508 + 666.0000 1.832765 0.7522531E-01 -60.48236 2319.509 + 668.0000 1.854995 0.7437853E-01 -60.48199 2312.552 + 670.0000 1.876976 0.7354498E-01 -60.48162 2305.636 + 672.0000 1.898712 0.7272440E-01 -60.48124 2298.762 + 674.0000 1.920205 0.7191652E-01 -60.48086 2291.928 + 676.0000 1.941461 0.7112110E-01 -60.48049 2285.135 + 678.0000 1.962482 0.7033788E-01 -60.48011 2278.381 + 680.0000 1.983272 0.6956663E-01 -60.47973 2271.668 + 682.0000 2.003834 0.6880710E-01 -60.47935 2264.994 + 684.0000 2.024173 0.6805907E-01 -60.47897 2258.358 + 686.0000 2.044291 0.6732231E-01 -60.47858 2251.762 + 688.0000 2.064193 0.6659661E-01 -60.47820 2245.204 + 690.0000 2.083880 0.6588175E-01 -60.47781 2238.683 + 692.0000 2.103356 0.6517751E-01 -60.47743 2232.201 + 694.0000 2.122625 0.6448371E-01 -60.47704 2225.755 + 696.0000 2.141689 0.6380014E-01 -60.47665 2219.347 + 698.0000 2.160551 0.6312659E-01 -60.47626 2212.976 + 700.0000 2.179215 0.6246290E-01 -60.47587 2206.640 diff --git a/tests/SetEpsTestCasesFromScratch/f90/extend3.f90 b/tests/SetEpsTestCasesFromScratch/f90/extend3.f90 new file mode 100644 index 0000000000000000000000000000000000000000..1ebde11e03092161284a5e2d949845e50271a1c0 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/extend3.f90 @@ -0,0 +1,22 @@ +subroutine extend3(array, extension) + +! increase the size of a allocatable array by extension. +! if the array is not allocated, create one of size extension. + + implicit none + + double precision, intent(in out), allocatable :: array(:, :) + integer, intent(in) :: extension + + integer :: old_size_1, old_size_2 + double precision, allocatable :: tmp_array(:, :) + + old_size_1 = size(array, 1) + old_size_2 = size(array, 2) + allocate(tmp_array(old_size_1, old_size_2 + extension)) + tmp_array(1:old_size_1, 1:old_size_2) = array + deallocate(array) + call move_alloc(tmp_array, array) + + return +end subroutine extend3 diff --git a/tests/SetEpsTestCasesFromScratch/f90/fint1.f90 b/tests/SetEpsTestCasesFromScratch/f90/fint1.f90 new file mode 100644 index 0000000000000000000000000000000000000000..62361a0b908e865c00fe498385a7e814d12e9941 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/fint1.f90 @@ -0,0 +1,64 @@ +module fint1_mod +contains +double precision function fint1(u, eps, thick, layers, nper, eps_size) + +! ****************************************************************** +! * * +! * integration over the azimutal angle from 0.0 to pi * +! * * +! * eps: epsilon * +! * thick: thickness * +! * layers: number of layers * +! * nper: number of periodic layers * +! * eps_size: number of layers * +! * * +! ****************************************************************** + + use surlos_mod + use param_mod + + implicit none + + double precision, intent(in) :: u + double complex, intent(in) :: eps(eps_size) + double precision, intent(in) :: thick(eps_size) + integer, intent(in) :: layers, nper, eps_size + + double precision :: den, dif, e, pi, rom, rop, sum, t, u2 + double precision :: usurlo + + data pi / 3.141592653589793238d0 / + +! write (*,*) 'fint1:' +! write (*,*) 'thick: ', size(thick) +! write (*,*) 'eps: ', size(eps) + + if (u == 0.0d0) then + fint1 = 0.0d0 + return + endif + e = tanpsi * u + u2 = u**2 + rom = (1.0d0 - e)**2 + u2 + rop = (1.0d0 + e)**2 + u2 + sum = rop + rom + rom = dsqrt(rom) + rop = dsqrt(rop) + dif = rop - rom + den = dsqrt((2.0d0 - dif) * (2.0d0 + dif)) * rop * rom + fint1 = pi * u2 * (4 * sum - dif**2 * (sum - rop * rom)) / den**3 + if (rational) then + return + endif + if (user) then + fint1 = fint1 * usurlo(ru * u, wn) + else + fint1 = fint1 * surlos(ru * u, eps, thick, layers, nper) + if (dlimf > 0.0d0) then + t = ru * u * dlimf + fint1 = fint1 * (1.d0 + t * dlog(t / (t + 0.26d0)))**2 / (1.d0 + 1.40d0 * t) + endif + endif + return +end function fint1 +end module fint1_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/fint2.f90 b/tests/SetEpsTestCasesFromScratch/f90/fint2.f90 new file mode 100644 index 0000000000000000000000000000000000000000..12baaf6b79c240f12c579dcc0d27e59add6e7dae --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/fint2.f90 @@ -0,0 +1,66 @@ +module fint2_mod +contains +double precision function fint2(u, eps, thick, layers, nper, eps_size) + +! ****************************************************************** +! * * +! * integration over the azimutal angle from 0.0 to phi < pi * +! * * +! * eps: epsilon * +! * thick: thickness * +! * layers: number of layers * +! * nper: number of periodic layers * +! * eps_size: number of layers * +! * * +! ****************************************************************** + + use surlos_mod + use param_mod + + implicit none + + double precision, intent(in) :: u + double complex, intent(in) :: eps(eps_size) + double precision, intent(in) :: thick(eps_size) + integer, intent(in) :: layers, nper, eps_size + + double precision :: arg, b2, c, phi, t, x + double precision :: phint, usurlo + +! write (*,*) 'fint2:' +! write (*,*) 'thick: ', size(thick) +! write (*,*) 'eps: ', size(eps) + + if (u == 0.0d0) then + fint2 = 0.0d0 + return + endif + b2 = bcoef**2 + c = (1.0d0 - elleps) * (cospsi * u)**2 + (bcoef - ccoef) * (bcoef + ccoef) + if (dabs(acoef * c) > 1.0d-03 * b2) then + x = (bcoef - dsqrt(b2 - acoef * c)) / acoef + else + x = acoef * c / b2 + x = 0.5d0 * c * (1.d0 + 0.25d0 * x * (1.d0 + 0.5d0 * x * (1.d0 + 0.625d0 * x))) / bcoef + endif + arg = x / u + if (dabs(arg) > 1.0d0) then + arg = dsign(1.0d0, arg) + endif + phi = dacos(arg) + fint2 = phint(phi, tanpsi, u) + if (rational) then + return + endif + if (user) then + fint2 = fint2 * usurlo(ru * u, wn) + else + fint2 = fint2 * surlos(ru * u, eps, thick, layers, nper) + if (dlimf > 0.0d0) then + t = ru * u * dlimf + fint2 = fint2 * (1.d0 + t * dlog(t / (t + 0.26d0)))**2 / (1.d0 + 1.40d0 * t) + endif + endif + return +end function fint2 +end module fint2_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/fint3.f90 b/tests/SetEpsTestCasesFromScratch/f90/fint3.f90 new file mode 100644 index 0000000000000000000000000000000000000000..7135e6f590a0acc464a660e71f853402881ef047 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/fint3.f90 @@ -0,0 +1,59 @@ +module fint3_mod +contains +double precision function fint3(u, eps, thick, layers, nper, eps_size) + +! ****************************************************************** +! * * +! * integration over the azimutal angle from phi1 > 0 to phi2 < pi * +! * * +! * eps: epsilon * +! * thick: thickness * +! * layers: number of layers * +! * nper: number of periodic layers * +! * eps_size: number of layers * +! * * +! ****************************************************************** + + use surlos_mod + use param_mod + + implicit none + + double precision, intent(in) :: u + double complex, intent(in) :: eps(eps_size) + double precision, intent(in) :: thick(eps_size) + integer, intent(in) :: layers, nper, eps_size + + double precision :: arg, phi1, phi2, rac, t + double precision :: phint, usurlo + +! write (*,*) 'fint3:' +! write (*,*) 'thick: ', size(thick) +! write (*,*) 'eps: ', size(eps) + + if (u == 0.0d0) then + fint3 = 0.0d0 + return + endif + rac = dsign(1.0d0, acoef) * cospsi * dsqrt((1.0d0 - elleps) * acoef * (um - u) * (um + u)) + arg = (bcoef - rac) / (u * acoef) + if (dabs(arg) > 1.0d0) arg = dsign(1.0d0, arg) + phi2 = dacos(arg) + fint3 = phint(phi2, tanpsi, u) + arg = (bcoef + rac) / (u * acoef) + if (dabs(arg) > 1.0d0) arg = dsign(1.0d0, arg) + phi1 = dacos(arg) + fint3 = fint3 - phint(phi1, tanpsi, u) + if (rational) return + if (user) then + fint3 = fint3 * usurlo(ru * u, wn) + else + fint3 = fint3 * surlos(ru * u, eps, thick, layers, nper) + if (dlimf > 0.0d0) then + t = ru * u * dlimf + fint3 = fint3 * (1.d0 + t * dlog(t / (t + 0.26d0)))**2 / (1.d0 + 1.40d0 * t) + endif + endif + return +end function fint3 +end module fint3_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/fun.f90 b/tests/SetEpsTestCasesFromScratch/f90/fun.f90 new file mode 100644 index 0000000000000000000000000000000000000000..b077a2adbb98b22677f9cc19855884d82b8333da --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/fun.f90 @@ -0,0 +1,23 @@ +double precision function fun(phi) + +! ****************************************************************** +! * * +! * integrand of the expression of the 1st order term in the * +! * expansion of the eels integral for a homogeneous target. * +! * * +! ****************************************************************** + + use param_mod + + implicit none + + double precision, intent(in) :: phi + + double precision :: sinphi + + sinphi = dsin(phi) + fun = dsqrt((1.0d0 - elleps + elleps * sinphi**2) * & + (1.0d0 - sinpsi * sinphi) * & + (1.0d0 + sinpsi * sinphi)) + return +end function fun diff --git a/tests/SetEpsTestCasesFromScratch/f90/get_commandline_options.f90 b/tests/SetEpsTestCasesFromScratch/f90/get_commandline_options.f90 new file mode 100644 index 0000000000000000000000000000000000000000..e9977be96a19336c3dbf6636d20485f8b56650ce --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/get_commandline_options.f90 @@ -0,0 +1,93 @@ +module get_commandline_options_mod +contains +subroutine get_commandline_options(eelsin_name, eelsou_name, bosin_name, bosou_name) + +! This routine defines the commandline options, parses the commandline and sets the +! filenames of I/O files. +! +! It uses sufr_getopt from libSUFR and has been derived from getopt_long_example.f90 +! +! KMS + + use sufr_getopt, only: getopt_t, getopt_long, longOption, optArg, getopt_long_help + + implicit none + + character (len = :), allocatable, intent(in out) :: eelsin_name, eelsou_name + character (len = :), allocatable, intent(in out) :: bosin_name, bosou_name + + include 'version.inc' + + character :: option + integer :: status + + ! Set up the longOpts struct to define the valid options: + ! short option, long option, argument (no = 0 / yes = 1), short description + type(getopt_t) :: longOpts(8) = & + [ & + getopt_t('v', 'version', 0, 'Print version'), & + getopt_t('V', 'version', 0, 'Print version'), & + getopt_t('h', 'help', 0, 'Print options'), & + getopt_t('d', 'dir', 1, 'Change I/O directory'), & + getopt_t('e', 'eelsin', 1, 'Name of EELS input file (eelsin)'), & + getopt_t('f', 'eelsou', 1, 'Name of EELS output file (optional)'), & + getopt_t('b', 'bosin', 1, 'Name of BOSON input file (bosin)'), & + getopt_t('c', 'bosou', 1, 'Name of BOSON output file (bosou)') & + ] + + do ! scan all the command-line parameters + + ! getopt_long() returns a single character" ">","!",".", or the short-option character (e.g. "a" for -a). + ! It also sets two 'global' variables through the SUFR_getopt module: + ! - longOption: the full option (e.g. "-a" or "--all") including the dashes + ! - optArg: the argument following the option (if required and present) + option = getopt_long(longOpts) + + ! Do different things depending on the option returned: + select case(option) + case('d') ! Change I/O directory + status = chdir(trim(optArg)) + if (status /= 0) then + write (*,*) 'WARNING: change directory failed!' + write (*,*) 'Directory tried: ', trim(optarg) + write (*,*) 'Error code (see: man chdir): ', status + write (*,*) 'Continuing in the start directory.' + write (*,*) '' + end if + case('e') + eelsin_name = optArg + case('f') + eelsou_name = optArg + case('b') + bosin_name = optArg + case('c') + bosou_name = optArg + case('v') + write (*,*) 'eels-boson version: ' // version + stop + case('V') + write (*,*) 'eels-boson version: ' // version + stop + case('h') + call getopt_long_help(longOpts) + stop + case('!') ! Unknown option (starting with "-" or "--") + write (*,*) 'WARNING: unknown option: ' // trim(optArg) + call getopt_long_help(longOpts) + stop + case('.') ! Parameter is not an option (i.e., it doesn't start with "-" or "--") + write (*,*) 'WARNING: parameter without option: ' // trim(optArg) + call getopt_long_help(longOpts) + stop + case default + write (*,*) 'WARNING: unknown option: ' // trim(longOption) + call getopt_long_help(longOpts) + stop + case('>') ! Last parameter. Exit case statement + exit + end select + end do + + return +end subroutine get_commandline_options +end module get_commandline_options_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/param.f90 b/tests/SetEpsTestCasesFromScratch/f90/param.f90 new file mode 100644 index 0000000000000000000000000000000000000000..a688b3cc98a5f11152df9c2432ada9863c71d093 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/param.f90 @@ -0,0 +1,5 @@ +module param_mod + double precision :: acoef, bcoef, ccoef, elleps, cospsi, sinpsi, tanpsi + double precision :: ru, um, dlimf, wn + logical :: user, rational +end module param_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/phint.f90 b/tests/SetEpsTestCasesFromScratch/f90/phint.f90 new file mode 100644 index 0000000000000000000000000000000000000000..57f15b5119872c2e2527c95ceac70168f8bfc8c8 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/phint.f90 @@ -0,0 +1,74 @@ +double precision function phint(phi, a, u) + +! ******************************************************************* +! * * +! * evaluate the integral from zero to phi of * +! * * +! * u 2 * +! * ( ----------------------------- ) dphi * +! * 2 2 * +! * (1 - a * u * cos(phi)) + u * +! * * +! * for 0 <= phi <= pi , u >= 0 and a >= 0 * +! * * +! * Reference: * +! * Ph. Lambin, J. P. Vigneron, and A. A. Lucas, * +! * Phys. Rev. B 32 (1985) 8203-8215. * +! * * +! ******************************************************************* + + implicit none + + double precision, intent(in) :: phi + double precision, intent(in) :: a + double precision, intent(in) :: u + + double precision :: ai, ar, bi, br, c, cpr, d, e, esr, pi, qr, ri, rm, root + double precision :: rp, rr, s, spr, tm, tp, u2, x, zeta, zetai, zetar, zr + + pi = 3.141592653589793238d0 + c = dcos(phi) + s = dsin(phi) + u2 = u**2 + e = a*u + if (u < 1.0d0 .and. e < 1.0d-02 * (1.0d0 + u2)) then + zr = 1.0d0 + u2 + esr = e / zr + phint = u2 / zr**2 * ((( (4.0d0 / 3.0d0) * (2.0d0 + c**2) * s * (5.0d0 - 3 * u2) * & + esr + (phi + c * s) * (5.0d0 - u2)) * esr + 4 * s) * esr + phi) + else + rm = dsqrt((1.0d0 - e)**2 + u2) + tm = 0.5d0 * datan2(u, 1.0d0 - e) + rp = dsqrt((1.0d0 + e)**2 + u2) + tp = 0.5d0 * datan2(u, 1.0d0 + e) + root = dsqrt(rm * rp) + cpr = dcos(tm + tp) + spr = dsin(tm + tp) + x = 0.0d0 ! ensure initialization + if (c >= 0.0d0) then + x = s / (1.0d0 + c) + elseif (dabs(s) > 1.0d-07) then + x = (1.0d0 - c) / s + endif + if ((c >= 0.0d0) .or. (dabs(s) > 1.0d-07)) then + zeta = dsqrt(rm / rp) + zetar = -zeta * dsin(tm - tp) + zetai = zeta * dcos(tm - tp) + br = 0.5d0 * dlog(((zetar + x)**2 + zetai**2) / ((zetar - x)**2 + zetai**2)) + bi = datan2(zetai, zetar + x) - datan2(zetai, zetar - x) + rr = -(br * spr - bi * cpr) / root + ri = -(bi * spr + br * cpr) / root + d = e * s / ((1.0d0 - e * c)**2 + u2) + ar = d * (1.0d0 - e * c) - rr + u * ri + ai = -d * u - ri - u * rr + else + rr = -pi / root * cpr + ri = pi / root * spr + ar = -rr + u * ri + ai = -ri - u * rr + endif + qr = (ar * (cpr - spr) * (cpr + spr) + 2 * ai * cpr * spr) / (rm * rp) + phint = 0.5d0 * (ri / u - qr) + endif + return +end function phint diff --git a/tests/SetEpsTestCasesFromScratch/f90/qrat.f90 b/tests/SetEpsTestCasesFromScratch/f90/qrat.f90 new file mode 100644 index 0000000000000000000000000000000000000000..7bbb966b70cc9d60ff1185d0fec5f9c293d5e29c --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/qrat.f90 @@ -0,0 +1,14 @@ +double precision function qrat(x, alpha, beta, c1, c2) + +! 1 + x * (beta + c1 * x) +! qrat = --------------------------------------------- +! (1 + x * (beta + c2 * x)) * (1 + alpha * x)^2 + + implicit none + + double precision, intent(in) :: x, alpha, beta, c1, c2 + + qrat = (1.0d0 + x * (beta + c1 * x)) / ((1.0d0 + x * (beta + c2 * x)) * (1.0d0 + alpha * x)**2) + + return +end function qrat diff --git a/tests/SetEpsTestCasesFromScratch/f90/quanc8.f90 b/tests/SetEpsTestCasesFromScratch/f90/quanc8.f90 new file mode 100644 index 0000000000000000000000000000000000000000..86ee4f36463a6f8150bb5f52b50528ca849879f5 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/quanc8.f90 @@ -0,0 +1,240 @@ +module quanc8_mod +contains +subroutine quanc8(fun, a, b, abserr, relerr, result, errest, nofun, flag, eps, thick, layers, nper) + +! estimate the integral of fun(x) from a to b +! to a user provided tolerance. +! an automatic adaptive routine based on +! the 8-panel newton-cotes rule (g. forsythe et al, 1977, p. 92) +! +! input .. +! +! fun the name of the integrand function subprogram fun(x). +! a the lower limit of integration. +! b the upper limit of integration.(b may be less than a.) +! relerr a relative error tolerance. (should be non-negative) +! abserr an absolute error tolerance. (should be non-negative) +! +! output .. +! +! result an approximation to the integral hopefully satisfying the +! least stringent of the two error tolerances. +! errest an estimate of the magnitude of the actual error. +! nofun the number of function values used in calculation of result. +! flag a reliability indicator. if flag is zero, then result +! probably satisfies the error tolerance. if flag is +! xxx.yyy , then xxx = the number of intervals which have +! not converged and 0.yyy = the fraction of the interval +! left to do when the limit on nofun was approached. +! +! additional input .. +! +! eps epsilon +! thick thickness +! layers number of layers +! nper number of periodic layers + + implicit none + + double precision :: fun + double precision, intent(in) :: a + double precision, intent(in) :: b + double precision, intent(in out) :: abserr + double precision, intent(in) :: relerr + double precision, intent(out) :: result + double precision, intent(out) :: errest + integer, intent(out) :: nofun + double precision, intent(out) :: flag + + external fun + + double precision, intent(in) :: thick(:) + double complex, intent(in) :: eps(:) + integer, intent(in) :: layers, nper + + double precision :: w0, w1, w2, w3, w4, area, x0, f0, stone, step, cor11, temp + double precision :: qprev, qnow, qdiff, qleft, esterr, tolerr + double precision :: qright(31), f(16), x(16), fsave(8, 30), xsave(8, 30) + double precision :: dabs, dmax1 + + integer :: levmin, levmax, levout, nomax, nofin, lev, nim, i, j + +! *** stage 1 *** general initialization +! set constants. + +! write (*,*) 'quanc8:' +! write (*,*) 'thick: ', size(thick) +! write (*,*) 'eps: ', size(eps) + + levmin = 1 + levmax = 30 + levout = 6 + nomax = 5000 + nofin = nomax - 8 * (levmax - levout + 2**(levout + 1)) + +! trouble when nofun reaches nofin + + w0 = 3956.0d0 / 14175.0d0 + w1 = 23552.0d0 / 14175.0d0 + w2 = -3712.0d0 / 14175.0d0 + w3 = 41984.0d0 / 14175.0d0 + w4 = -18160.0d0 / 14175.0d0 + +! initialize running sums to zero. + + flag = 0.0d0 + result = 0.0d0 + cor11 = 0.0d0 + errest = 0.0d0 + area = 0.0d0 + nofun = 0 + if (a == b) return + +! *** stage 2 *** initialization for first interval + + lev = 0 + nim = 1 + x0 = a + x(16) = b + qprev = 0.0d0 + f0 = fun(x0, eps, thick, layers, nper, size(eps)) + stone = (b - a) / 16 + x(8) = (x0 + x(16)) / 2 + x(4) = (x0 + x(8)) / 2 + x(12) = (x(8) + x(16)) / 2 + x(2) = (x0 + x(4)) / 2 + x(6) = (x(4) + x(8)) / 2 + x(10) = (x(8) + x(12)) / 2 + x(14) = (x(12) + x(16)) / 2 +!*** !$OMP PARALLEL shared (f, x, eps, thick, layers, nper) private (j) +!*** !$OMP DO + do j = 2, 16, 2 + f(j) = fun(x(j), eps, thick, layers, nper, size(eps)) + enddo +!*** !$OMP END DO NOWAIT +!*** !$OMP END PARALLEL + nofun = 9 + + do + +! *** stage 3 *** central calculation +! requires qprev, x0, x2, x4, ..., x16, f0, f2, f4, ..., f16. +! calculates x1, x3, ...x15, f1, f3, ...f15, qleft, qright, qnow, qdiff, area. + + x(1) = (x0 + x(2)) / 2 + f(1) = fun(x(1), eps, thick, layers, nper, size(eps)) + do j = 3, 15, 2 + x(j) = (x(j - 1) + x(j + 1)) / 2 + enddo +!*** !$OMP PARALLEL shared (f, x, eps, thick, layers, nper) private (j) +!*** !$OMP DO SCHEDULE(auto) + do j = 3, 15, 2 + f(j) = fun(x(j), eps, thick, layers, nper, size(eps)) + enddo +!*** !$OMP END DO NOWAIT +!*** !$OMP END PARALLEL + nofun = nofun + 8 + step = (x(16) - x0) / 16.0d0 + qleft = (w0 * (f0 + f(8)) + w1 * (f(1) + f(7)) + w2 * (f(2) + f(6)) & + + w3 * (f(3) + f(5)) + w4 * f(4)) * step + qright(lev + 1) = (w0 * (f(8) + f(16)) + w1 * (f(9) + f(15)) + w2 * (f(10) + f(14)) & + + w3 * (f(11) + f(13)) + w4 * f(12)) * step + qnow = qleft + qright(lev + 1) + qdiff = qnow - qprev + area = area + qdiff + +! *** stage 4 *** interval convergence test + + esterr = dabs(qdiff) / 1023 + tolerr = dmax1(abserr, relerr * dabs(area)) * (step / stone) + + if (lev >= levmin) then + if (lev >= levmax) then +! current level is levmax. + flag = flag + 1.0d0 + else + if (nofun > nofin) then +! *** stage 6 *** trouble section +! number of function values is about to exceed limit. + nofin = 2 * nofin + levmax = levout + flag = flag + (b - x0) / (b - a) + else + if (esterr > tolerr) then +! *** stage 5 *** no convergence +! locate next interval. + nim = 2 * nim + lev = lev + 1 +! store right hand elements for future use. + do i = 1, 8 + fsave(i, lev) = f(i + 8) + xsave(i, lev) = x(i + 8) + enddo +! assemble left hand elements for immediate use. + qprev = qleft + do i = 1, 8 + f(18 - 2 * i) = f(9 - i) + x(18 - 2 * i) = x(9 - i) + enddo + cycle + endif + endif + endif + +! *** stage 7 *** interval converged +! add contributions into running sums. + result = result + qnow + errest = errest + esterr + cor11 = cor11 + qdiff / 1023 +! locate next interval. + do while (nim /= 2 * (nim / 2)) + nim = nim / 2 + lev = lev - 1 + enddo + nim = nim + 1 + + if (lev <= 0) exit + +! assemble elements required for the next interval. + qprev = qright(lev) + x0 = x(16) + f0 = f(16) + do i = 1, 8 + f(2*i) = fsave(i, lev) + x(2*i) = xsave(i, lev) + enddo + cycle + else +! *** stage 5 *** no convergence +! locate next interval. + nim = 2 * nim + lev = lev + 1 +! store right hand elements for future use. + do i = 1, 8 + fsave(i, lev) = f(i + 8) + xsave(i, lev) = x(i + 8) + enddo +! assemble left hand elements for immediate use. + qprev = qleft + do i = 1, 8 + f(18 - 2 * i) = f(9 - i) + x(18 - 2 * i) = x(9 - i) + enddo + endif + + enddo + + ! *** stage 8 *** finalize and return + result = result + cor11 + +! make sure errest not less than roundoff level. + if (errest /= 0.0d0) then + temp = dabs(result) + errest + do while (temp == dabs(result)) + errest = 2 * errest + temp = dabs(result) + errest + enddo + endif + return +end subroutine quanc8 +end module quanc8_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/queels.f90 b/tests/SetEpsTestCasesFromScratch/f90/queels.f90 new file mode 100644 index 0000000000000000000000000000000000000000..d2128df921e17e578004094b2218052d79c1cbe1 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/queels.f90 @@ -0,0 +1,92 @@ +module queels_mod +contains +subroutine queels(x, f, aerr, rerr, facru, eps, thick, layers, nper) + +! ****************************************************************** +! * * +! * perform q-space integration for computing the eels spectrum of * +! * a isotropic target using polar coordinates. * +! * * +! * x is the dimensionless energy loss hbar*omega/(2*e0*phia) * +! * aerr and rerr are the desired absolute and relative accuracies * +! * facru*x is the units of wavevectors omega/v_perpendicular * +! * f is the q-integral multiplied by (2/pi)**2 * +! * * +! * eps: epsilon * +! * thick: thickness * +! * layers: number of layers * +! * nper: number of periodic layers * +! * * +! ****************************************************************** + + use quanc8_mod + use fint1_mod + use fint2_mod + use fint3_mod + use param_mod + + implicit none + + double precision, intent(in) :: x + double precision, intent(out) :: f + double precision, intent(in out) :: aerr + double precision, intent(in out) :: rerr + double precision, intent(in) :: facru + double complex, intent(in) :: eps(:) + double precision, intent(in) :: thick(:) + integer, intent(in) :: layers, nper + + double precision :: error, flag + double precision :: u1, u2, ut, y + integer :: ie, nofu + dimension error(3), flag(3) + +! write (*,*) 'queels:' +! write (*,*) 'eps: ', size(eps) +! write (*,*) 'thick: ', size(thick) + + f = 0.0d0 + if (x <= 0.0d0) then + return + endif + ru = facru * x + ccoef = cospsi**2 / x + ut = ccoef - bcoef + u1 = dabs(ut) + u2 = ccoef + bcoef + if (ut > 0.0d0) then + call quanc8(fint1, 0.0d0, u1, aerr, rerr, y, error(1), nofu, flag(1), eps, thick, layers, nper) + f = y + else + flag(1) = 0.0d0 + endif + if (u2 > u1) then + call quanc8(fint2, u1, u2, aerr, rerr, y, error(2), nofu, flag(2), eps, thick, layers, nper) + f = f + y + else + flag(2) = 0.0d0 + endif + if (dabs(acoef) > x * (1.0d0 - elleps) * bcoef) then + um = dsqrt(ccoef / x / (1.0d0 - elleps) + bcoef**2 / acoef) + if (um > u2) then + call quanc8(fint3, u2, um, aerr, rerr, y, error(3), nofu, flag(3), eps, thick, layers, nper) + f = f + y + endif + if (um < u1) then + call quanc8(fint3, um, u1, aerr, rerr, y, error(3), nofu, flag(3), eps, thick, layers, nper) + f = f - y + endif + else + flag(3) = 0.0d0 + endif + do ie = 1, 3 + if (flag(ie) == 0.0d0) cycle + write(*,*) ' +++ flag(', ie, ') =', flag(ie), ', error =', error(ie), ' +++' + if (flag(ie) - aint(flag(ie)) > 0.5d-02) then + stop '*** execution aborted ***' + endif + enddo + f = (2 / 3.141592653589793238d0)**2 * f + return +end subroutine queels +end module queels_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/seteps.f90 b/tests/SetEpsTestCasesFromScratch/f90/seteps.f90 new file mode 100644 index 0000000000000000000000000000000000000000..323986b84876e1b4ef80715e5beea58adb0a6a56 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/seteps.f90 @@ -0,0 +1,103 @@ +module seteps_mod +contains +subroutine seteps(neps, nos, osc, epsinf, wn, name, eps, layers, mode) + +! ****************************************************************** +! * * +! * set up long-wavelength dielectric functions of the layers for * +! * the present frequency wn (in cm**-1) * +! * * +! * neps: number of epsilon values * +! * nos: number of oscillators * +! * osc: oscillator parameters, wTO, Q, lambda * +! * epsinf: epsilon at infinite frequency * +! * wn: frequency * +! * name: layer names * +! * eps: epsilon * +! * layers: number of layers * +! * mode: 'kurosawa' for kurosawa model * +! * * +! ****************************************************************** + + implicit none + integer, intent(in) :: neps + integer, intent(in) :: nos(:) + double precision, intent(in) :: osc(:, :), epsinf(:), wn + character (len=10), intent(in) :: name(:) + double complex, intent(in out) :: eps(:) + integer, intent(in) :: layers + character (len=10), intent(in) :: mode + + double precision :: x + double complex :: deno, nomi + integer :: j, k, l, m + +! write (*,*) 'seteps:' +! write (*,*) 'nos: ', size(nos) +! write (*,*) 'osc: ', size(osc) +! write (*,*) 'epsinf: ', size(epsinf) +! write (*,*) 'name: ', size(name) +! write (*,*) 'eps: ', size(eps) +! write (*,*) 'thick: ', size(thick) + + j = 0 + do l = 1, neps + m = nos(l)/2 + nomi = dcmplx(1.0d0, 0.0d0) + deno = dcmplx(1.0d0, 0.0d0) + if (mode == 'kurosawa') then + do k = 1, m + j = j + 1 +! since osc(1,*) and wn are real, the following should be equivalent +! Check required. +! Furthermore the first term can be rewritten as +! (osc(1, j + m) - wn) * (osc(1, j + m) + wn) +! nomi = nomi * cmplx(osc(1, j + m)**2 - wn**2, - wn * osc(3, j + m)) + nomi = nomi * (osc(1, j + m)**2 - wn**2 - dcmplx(0.0d0, wn * osc(3, j + m))) + deno = deno * (osc(1, j )**2 - wn**2 - dcmplx(0.0d0, wn * osc(3, j ))) + enddo + eps(l) = epsinf(l) * nomi / deno + elseif (name(l) == 'metal') then + j = j + 1 +! suggestion for replacement +! eps(l) = -osc(1, j)**2 / cmplx(wn**2, wn * osc(3, j)) + eps(l) = -osc(1, j)**2 / ( wn**2 + dcmplx(0.0d0, wn * osc(3, j)) ) + else + eps(l) = epsinf(l) +! The original version had this additional loop. It seems, it has been removed +! because all cases of nos(l) > 1 are treated in case 1 above + do k = 1, nos(l) + j = j + 1 + if (osc(1, j) > 1d-3) then ! prevent division by zero + x = wn / osc(1, j) + else + x = wn * 1d3 + endif + deno = x * dcmplx(x, osc(3, j)) + if (osc(2, j) >= 0.0d0) then + deno = 1.0d0 - deno + endif + if (cdabs(deno) == 0.0d0) then ! replace 0 by machine epsilon + ! if deno is always > 0 then this would do it: + ! deno = cdmax(deno, epsilon(1.0d0) / 2) + deno = epsilon(1.0d0) / 2 + endif + eps(l) = eps(l) + osc(2, j) / deno + enddo + endif + enddo + if (neps == layers + 1) then +! the substrate is a anisotropic uniaxial material + eps(layers) = cdsqrt(eps(layers) * eps(layers + 1)) + if (dimag(eps(layers)) < 0.0d0) then + eps(layers) = -eps(layers) + endif + endif + +! **** log modification start + write (99, '(30g15.7)') wn, (eps(j), j = 1, neps) +! **** log modification end + + return +end subroutine seteps +end module seteps_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/surlos.f90 b/tests/SetEpsTestCasesFromScratch/f90/surlos.f90 new file mode 100644 index 0000000000000000000000000000000000000000..66552e98a2e78ab9c2a4cd0b42c08d7a457c9a4e --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/surlos.f90 @@ -0,0 +1,174 @@ +module surlos_mod +contains + +logical function zero(z) + double complex, intent(in) :: z + zero = (dble(z) == 0.0d0) .and. (dimag(z) == 0.0d0) +end function zero + +double precision function surlos(dk, eps, thick, layers, nper) + +! ****************************************************************** +! * * +! * eels surface loss function for an arbitrary multilayered target* +! * * +! * dk: * +! * eps: epsilon * +! * thick: thickness * +! * layers: number of layers * +! * nper: number of periodic layers * +! * * +! ****************************************************************** + + implicit none + + double precision, intent(in) :: dk + double complex, intent(in) :: eps(:) + double precision, intent(in) :: thick(:) + integer, intent(in) :: layers, nper + + integer :: lmax + logical :: static, skip + integer :: lstart, n + double precision, allocatable :: arg(:) + double precision :: argmin, argmax, cn, cnm1, epsmac, sn, snm1, tn + double complex :: a, b, csi, pnm2, pnm1, pn, pp, qnm2, qnm1, qn, qp + + epsmac = epsilon(1.0d0) + argmin = dsqrt(2 * epsmac) + epsmac = epsmac / 2 + argmax = 0.5d0 * dlog(2 / epsmac) + +! write (*,*) 'surlos:' +! write (*,*) 'thick: ', size(thick) +! write (*,*) 'eps: ', size(eps) + + lmax = size(eps) + allocate (arg(lmax)) + lstart = layers - nper + 1 + static = .true. + skip = .false. + + do n = 1, layers + arg(n) = dk * thick(n) + if (arg(n) > argmax .or. zero(eps(n))) then + csi = eps(n) + skip = .true. + exit + endif + static = .not. (n >= lstart .and. arg(n) > argmin) + enddo + + if (.not. skip) then +! *** periodic continued fraction, period = nper + + do ! Do loop is only for the option to skip out of it. + if (nper <= 1) then + csi = eps(layers) + else + if (static) then + pn = 0.0d0 + qn = 0.0d0 + do n = lstart, layers + pn = pn + thick(n) * eps(n) + qn = qn + thick(n) / eps(n) + enddo + if (zero(qn)) then + n = lstart + if (n <= 1) then + surlos = 0.0d0 + return + endif + n = n - 1 + csi = eps(n) / dtanh(arg(n)) + exit + endif ! zero(qn) + csi = cdsqrt(pn / qn) + if ((dimag(csi) < 0.0d0) .and. (dble(qn) < 0.0d0)) then + csi = -csi + endif + else ! static + cn = dcosh(arg(lstart)) + sn = dsinh(arg(lstart)) + pnm1 = 1.0d0 + pn = cn + pp = eps(lstart) * sn + qnm1 = 0.0d0 + qn = sn / eps(lstart) + qp = pn + do n = lstart + 1, layers + cnm1 = cn + snm1 = sn + cn = dcosh(arg(n)) + sn = dsinh(arg(n)) + a = eps(n) * sn + pp = cn * pp + a * pn + qp = cn * qp + a * qn + b = (eps(n - 1) / eps(n)) * (sn / snm1) + a = cnm1 * b + cn + pnm2 = pnm1 + pnm1 = pn + qnm2 = qnm1 + qnm1 = qn + pn = a * pnm1 - b * pnm2 + qn = a * qnm1 - b * qnm2 + enddo + if (zero(qn)) then + a = qp - pn + if (zero(a)) then + n = lstart + if (n <= 1) then + surlos = 0.0d0 + return + endif + n = n - 1 + csi = eps(n) / dtanh(arg(n)) + exit + endif + csi = pp / a + else + a = 0.5d0 * (pn - qp) / qn + b = cdsqrt(a**2 + pp / qn) + pn = a - pn / qn + if (cdabs(pn + b) > cdabs(pn - b)) then + b = -b + endif + csi = a + b + endif ! zero(qn) + endif ! static +! *** small-dk limit of the periodic tail + + endif ! nper <= 1 + + n = lstart + + exit + enddo + endif ! .not. skip + +! *** backward algorithm + + do + n = n - 1 + if (n <= 0) then + a = csi + 1.0d0 + if (zero(a)) then + surlos = 2 / epsmac + else + surlos = dimag(-2 / a) + endif + return + endif ! n <= 0 + + if (arg(n) /= 0.0d0) then + tn = dtanh(arg(n)) + b = eps(n) + csi * tn + if (zero(b)) then + surlos = 0.0d0 + return + endif + csi = eps(n) * (csi + tn * eps(n)) / b + endif + enddo +end function surlos +end module surlos_mod diff --git a/tests/SetEpsTestCasesFromScratch/f90/usurlo.f90 b/tests/SetEpsTestCasesFromScratch/f90/usurlo.f90 new file mode 100644 index 0000000000000000000000000000000000000000..222946c95f302c8015f157572135fbf0d44c1174 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/usurlo.f90 @@ -0,0 +1,23 @@ +double precision function usurlo(dq, wn) + +! ****************************************************************** +! * * +! * user-supplied dielectric surface loss function aimag(g(dq, wn)) * +! * input arguments : * +! * dq : modulus of the two-dimensional surface wave vector * +! * (angstroem**-1) * +! * wn : frequency (cm**-1) * +! * * +! ****************************************************************** + + implicit none + + double precision, intent(in out) :: dq + double precision, intent(in out) :: wn + + if ((wn .GT. 0) .AND. (dq .GT. 0)) then + write(*,*) 'hello, here is the user loss function' + endif + usurlo = 1.0d0 + return +end function usurlo diff --git a/tests/SetEpsTestCasesFromScratch/f90/version.inc b/tests/SetEpsTestCasesFromScratch/f90/version.inc new file mode 100644 index 0000000000000000000000000000000000000000..453b451b102cb59869bcbda41f7d309b7459ada2 --- /dev/null +++ b/tests/SetEpsTestCasesFromScratch/f90/version.inc @@ -0,0 +1,2 @@ +! version of eels-boson + character (len = *), parameter :: version = '1.0.0'