From b9823c0c90e4c4d98ac779943bb7696fb50ed8ff Mon Sep 17 00:00:00 2001
From: kamischi <karl-michael.schindler@web.de>
Date: Wed, 8 May 2024 17:25:52 +0200
Subject: [PATCH] refactor common / mulayr /
remove common / mulayr / by
- replacing the calculation of epsmac by call of function epsilon
- move calculation of argmin and argmax to the function surlos. They are used there only.
---
source/f90/doeels.f90 | 11 ++---------
source/f90/seteps.f90 | 8 +++-----
source/f90/surlos.f90 | 5 ++++-
3 files changed, 9 insertions(+), 15 deletions(-)
diff --git a/source/f90/doeels.f90 b/source/f90/doeels.f90
index 40f39fb..d953363 100644
--- a/source/f90/doeels.f90
+++ b/source/f90/doeels.f90
@@ -55,8 +55,8 @@ subroutine doeels (e0, theta, phia, phib, wmin, wmax, dw, comment, comment_size,
logical :: rational, user
integer :: i, iw, neps, nofu, nout, nw, lmax
- double precision :: a, acoef, aerr, alpha, argmin, argmax, b, bcoef, beta, &
- c1, c2, ccoef, cospsi, dlimf, dx, elleps, ener, epsmac, facru, f, f0, &
+ double precision :: a, acoef, aerr, alpha, b, bcoef, beta, &
+ c1, c2, ccoef, cospsi, dlimf, dx, elleps, ener, facru, f, f0, &
f1, fpic, fun, pi, prefac, psia, psii, qrat, rerr, ru, sinpsi, t, &
tanpsi, table, um, widt, wn, wpic, x, xmin, xmax, z, z1, z2
double complex, allocatable :: eps(:)
@@ -66,7 +66,6 @@ subroutine doeels (e0, theta, phia, phib, wmin, wmax, dw, comment, comment_size,
common / param / acoef, bcoef, ccoef, elleps, cospsi, sinpsi, tanpsi, &
ru, um, dlimf, wn, user, rational
- common / mulayr / argmin, argmax, epsmac
data aerr / 0.0d0 /, rerr / 1.0d-06 /, f / 0.0d0 /, f1 / 0.0d0 /
@@ -86,12 +85,6 @@ subroutine doeels (e0, theta, phia, phib, wmin, wmax, dw, comment, comment_size,
! *** epsmac + 1.0 = epsmac , cosh(argmin) = 1.0 , tanh(argmax) = 1.0
pi = 4 * datan(1.0d0)
- epsmac = 1.0d0
- do while (1.0d0 + epsmac > 1.0d0)
- epsmac = epsmac / 2
- enddo
- argmin = dsqrt(2 * epsmac)
- argmax = 0.5d0 * dlog(2 / epsmac)
dlimf = 0.0d0
rational = .false.
diff --git a/source/f90/seteps.f90 b/source/f90/seteps.f90
index e9604f4..2900bcd 100644
--- a/source/f90/seteps.f90
+++ b/source/f90/seteps.f90
@@ -31,12 +31,10 @@ subroutine seteps(neps, nos, osc, epsinf, wn, name, eps, layers, mode)
double complex, intent(in out) :: eps(:)
integer, intent(in) :: layers
- double precision :: argmin, argmax, epsmac, x
+ double precision :: x
double complex :: deno, nomi
integer :: j, k, l, m
- common / mulayr / argmin, argmax, epsmac
-
! write (*,*) 'seteps:'
! write (*,*) 'nos: ', size(nos)
! write (*,*) 'osc: ', size(osc)
@@ -78,8 +76,8 @@ subroutine seteps(neps, nos, osc, epsinf, wn, name, eps, layers, mode)
if (osc(2, j) >= 0.0d0) then
deno = 1.0d0 - deno
endif
- if (cdabs(deno) == 0.0d0) then
- deno = epsmac
+ if (cdabs(deno) == 0.0d0) then ! replace 0 by machine epsilon
+ deno = epsilon(1.0d0) / 2
endif
eps(l) = eps(l) + osc(2, j) / deno
enddo
diff --git a/source/f90/surlos.f90 b/source/f90/surlos.f90
index 71bd2e9..ec92d72 100644
--- a/source/f90/surlos.f90
+++ b/source/f90/surlos.f90
@@ -28,7 +28,10 @@ double precision function surlos(dk, eps, thick, layers, nper)
double precision :: argmin, argmax, cn, cnm1, epsmac, sn, snm1, tn
double complex :: a, b, csi, pnm2, pnm1, pn, pp, qnm2, qnm1, qn, qp, z
- common / mulayr / argmin, argmax, epsmac
+ epsmac = epsilon(1.0d0)
+ argmin = dsqrt(2 * epsmac)
+ epsmac = epsmac / 2
+ argmax = 0.5d0 * dlog(2 / epsmac)
zero(z) = (dble(z) == 0.0d0) .and. (dimag(z) == 0.0d0)
--
GitLab