[FLASH-USERS] Compiling LaserSlab for FLASH 4.4 error in Diffuse_computeFluxLimiter.F90

[Klaus Weide]

On Fri, 30 Dec 2016, Marissa Adams wrote:

> mpixlf90_r -I/usr/local/hypre/2.8.0b-MPI-XL-No-Global-Partition/include -g
> > -O3 -qnohot -qrealsize=8 -qnosave -qfixed -c -qthreaded
> > -qsuffix=f=F90:cpp=F90 -qfree=f90 -WF,-DFLASH_3T -WF,-DMAXBLOCKS=1000
> > -WF,-DNXB=16 -WF,-DNYB=16 -WF,-DNZB=1 -WF,-DN_DIM=2
> > Diffuse_computeFluxLimiter.F90
> > "Diffuse_computeFluxLimiter.F90", line 200.37: 1513-041 (S) Arguments of
> > the wrong type were specified for the INTRINSIC procedure "min".
> > ** diffuse_computefluxlimiter   === End of Compilation 1 ===
> > 1501-511  Compilation failed for file Diffuse_computeFluxLimiter.F90.
> > make: *** [Diffuse_computeFluxLimiter.o] Error 1
> 
> 
> ​The line in question is:
> 
> dcoef = min(dcoefOld, fl/(maggrad + 1.0D-100))

Apply the attached patch (patch -p0 < r25554.diff in the top directory 
should do it).


> ​It seems to me that min does not like the mixture of two data types:
> variables *dcoefold*, *fl*, *maggrad* being reals, and a double precision
> number. I tried converting to and from double precision, and it still
> yields the same error but up upward in the file at line 162.40, which
> doesn't include a *min* at all. Gah!


It's a problem because different compilers handle the automatic 
propagation of Fortran default REAL to 8-byte REAL differently.
We could make the problem go away for IBM's XLF compiler by replacing
all 'D' exponents in real constants by 'E'. But then some other
compiler complains (for values outside the 4-byte IEEE single 
precision real range).

The approach taken in the patch should work in the maximum number
of compilers.

Klaus
-------------- next part --------------
Index: source/physics/Diffuse/DiffuseMain/Diffuse_computeFluxLimiter.F90
===================================================================
--- source/physics/Diffuse/DiffuseMain/Diffuse_computeFluxLimiter.F90	(revision 25553)
+++ source/physics/Diffuse/DiffuseMain/Diffuse_computeFluxLimiter.F90	(revision 25554)
@@ -86,6 +86,7 @@
   integer, intent(IN) :: blockID
   integer, intent(IN),OPTIONAL :: gcLayers
 
+  integer, parameter :: r8 = kind(1.0)
   integer :: lb, i, j, k
   integer :: blklim(2,MDIM), blklimgc(2,MDIM)
   integer :: ng
@@ -193,15 +194,15 @@
 
            select case(mode)
            case(FL_HARMONIC)
-              dcoef = 1.0 / (1.0/dcoefOld + maggrad/(fl + 1.0D-100))
+              dcoef = 1.0 / (1.0/dcoefOld + maggrad/(fl + 1.0e-100_r8))
               lambda3 = dcoef / dcoefOld
 
            case(FL_MINMAX)
-              dcoef = min(dcoefOld, fl/(maggrad + 1.0d-100))
+              dcoef = min(dcoefOld, fl/(maggrad + 1.0e-100_r8))
               lambda3 = dcoef / dcoefOld
 
            case(FL_LARSEN)
-              dcoef = 1.0 / sqrt((1.0/dcoefOld)**2 + (maggrad/(fl + 1.0D-100))**2)
+              dcoef = 1.0 / sqrt((1.0/dcoefOld)**2 + (maggrad/(fl + 1.0e-100_r8))**2)
               lambda3 = dcoef / dcoefOld
 
            case(FL_LEVPOM)
Index: source/physics/Diffuse/DiffuseMain/Diffuse_fluxLimiter.F90
===================================================================
--- source/physics/Diffuse/DiffuseMain/Diffuse_fluxLimiter.F90	(revision 25553)
+++ source/physics/Diffuse/DiffuseMain/Diffuse_fluxLimiter.F90	(revision 25554)
@@ -73,6 +73,7 @@
   integer, intent(IN) :: blkcnt
   integer, intent(IN) :: blklst(blkcnt)
 
+  integer, parameter :: r8 = kind(1.0)
   integer :: lb, i, j, k
   real, pointer :: blkPtr(:,:,:,:)
   integer :: blklim(2,MDIM), blklimgc(2,MDIM)
@@ -155,15 +156,15 @@
 
               select case(mode)
               case(FL_HARMONIC)
-                 dcoef = 1.0 / (1.0/dcoefOld + maggrad/(fl + 1.0D-100))
+                 dcoef = 1.0 / (1.0/dcoefOld + maggrad/(fl + 1.0e-100_r8))
                  lambda3 = dcoef / dcoefOld
 
               case(FL_MINMAX)
-                 dcoef = min(dcoefOld, fl/(maggrad + 1.0d-100))
+                 dcoef = min(dcoefOld, fl/(maggrad + 1.0e-100_r8))
                  lambda3 = dcoef / dcoefOld
 
               case(FL_LARSEN)
-                 dcoef = 1.0 / sqrt((1.0/dcoefOld)**2 + (maggrad/(fl + 1.0D-100))**2)
+                 dcoef = 1.0 / sqrt((1.0/dcoefOld)**2 + (maggrad/(fl + 1.0e-100_r8))**2)
                  lambda3 = dcoef / dcoefOld
 
               case(FL_LEVPOM)