diff --git a/src/FEQParse.F90 b/src/FEQParse.F90
index 7b76be2..c3564da 100644
--- a/src/FEQParse.F90
+++ b/src/FEQParse.F90
@@ -270,7 +270,7 @@ subroutine Tokenize(parser,tokenized,errorMsg)
 
         i = i + 1
 
-      elseif (parser % func % IsFunction(parser % inFixFormula(i:i))) then
+      elseif (IsFunction(parser % inFixFormula(i:i))) then
 
         parser % inFix % top_index = parser % inFix % top_index + 1
         parser % inFix % tokens(parser % inFix % top_index) % tokenString = ''
@@ -1673,7 +1673,7 @@ integer function Priority(parser,operatorString)
     class(EquationParser) :: parser
     character(*) :: operatorString
 
-    if (parser % func % IsFunction(operatorString)) then
+    if (IsFunction(operatorString)) then
 
       Priority = 4
 
@@ -1696,5 +1696,28 @@ integer function Priority(parser,operatorString)
     end if
 
   end function Priority
+  logical function IsFunction(eqChar)
+  character(*) :: eqChar
+
+  IsFunction = .false.
+  if (eqChar(1:1) == "\") then
+    IsFunction = .true.
+  end if
+
+end function IsFunction
+
+function FindLastFunctionIndex(eqChar) result(j)
+  character(*) :: eqChar
+  integer      :: i,j
+
+  do i = 1,len(eqChar)
+    if (eqChar(i:i) == "(") then
+      j = i - 2
+      exit
+    end if
+
+  end do
+
+end function FindLastFunctionIndex
 
 end module FEQParse
diff --git a/src/FEQParse_Functions.F90 b/src/FEQParse_Functions.F90
index 5fc72b7..db3ddf9 100644
--- a/src/FEQParse_Functions.F90
+++ b/src/FEQParse_Functions.F90
@@ -30,7 +30,7 @@ module FEQParse_Functions
 
   contains
     procedure,public :: Destruct => Destruct_FEQParse_FunctionHandler
-    procedure,public :: IsFunction
+    !procedure,public :: IsFunction
     generic,public :: f_of_x => f_of_x_sfp32, &
       f_of_x_r1fp32, &
       f_of_x_r2fp32, &
@@ -92,37 +92,6 @@ subroutine Destruct_FEQParse_FunctionHandler(functionhandler_obj)
 
   end subroutine Destruct_FEQParse_FunctionHandler
 
-  logical function IsFunction(functionhandler_obj,eqChar)
-    class(FEQParse_FunctionHandler) :: functionhandler_obj
-    character(*) :: eqChar
-    ! Local
-    integer :: i
-
-    IsFunction = .false.
-    do i = 1,functionhandler_obj % nFunctions
-
-      if (eqChar(1:1) == "\") then
-        IsFunction = .true.
-      end if
-
-    end do
-
-  end function IsFunction
-
-  function FindLastFunctionIndex(eqChar) result(j)
-    character(*) :: eqChar
-    integer      :: i,j
-
-    do i = 1,len(eqChar)
-      if (eqChar(i:i) == "(") then
-        j = i - 2
-        exit
-      end if
-
-    end do
-
-  end function FindLastFunctionIndex
-
   subroutine f_of_x_sfp32(functionhandler_obj,func,x,fx)
     !! Evaluates function for scalar fp32 input and output
     class(FEQParse_FunctionHandler) :: functionhandler_obj
diff --git a/test/abs_r3fp32.f90 b/test/abs_r3fp32.f90
index d5d9245..735ae35 100644
--- a/test/abs_r3fp32.f90
+++ b/test/abs_r3fp32.f90
@@ -27,14 +27,14 @@ integer function abs_r3fp32() result(r)
 
   x = 0.0_real32
   do k = 1,N
-  do j = 1,N
-    do i = 1,N
-      x(i,j,k,1) = -1.0_real32 + (2.0_real32)/real(N,real32)*real(i - 1,real32)
-      x(i,j,k,2) = -1.0_real32 + (2.0_real32)/real(N,real32)*real(j - 1,real32)
-      x(i,j,k,3) = -1.0_real32 + (2.0_real32)/real(N,real32)*real(k - 1,real32)
+    do j = 1,N
+      do i = 1,N
+        x(i,j,k,1) = -1.0_real32 + (2.0_real32)/real(N,real32)*real(i - 1,real32)
+        x(i,j,k,2) = -1.0_real32 + (2.0_real32)/real(N,real32)*real(j - 1,real32)
+        x(i,j,k,3) = -1.0_real32 + (2.0_real32)/real(N,real32)*real(k - 1,real32)
+      end do
     end do
   end do
-  end do
 
   do k = 1,N
     do j = 1,N
@@ -46,7 +46,7 @@ integer function abs_r3fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= 10.0_real32*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/abs_r3fp64.f90 b/test/abs_r3fp64.f90
index 75280f9..3f5b87c 100644
--- a/test/abs_r3fp64.f90
+++ b/test/abs_r3fp64.f90
@@ -46,7 +46,7 @@ integer function abs_r3fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= 10.0_real64*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/abs_r4fp32.f90 b/test/abs_r4fp32.f90
index 3655ad2..5e13751 100644
--- a/test/abs_r4fp32.f90
+++ b/test/abs_r4fp32.f90
@@ -51,9 +51,10 @@ integer function abs_r4fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <=  maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
+    print*,maxval(abs(feval - fexact))
     r = 1
   end if
 
diff --git a/test/abs_r4fp64.f90 b/test/abs_r4fp64.f90
index dff5c1e..5b6d93f 100644
--- a/test/abs_r4fp64.f90
+++ b/test/abs_r4fp64.f90
@@ -51,7 +51,7 @@ integer function abs_r4fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/acos_r3fp32.f90 b/test/acos_r3fp32.f90
index b06fdaf..67da284 100644
--- a/test/acos_r3fp32.f90
+++ b/test/acos_r3fp32.f90
@@ -46,7 +46,7 @@ integer function acos_r3fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/acos_r3fp64.f90 b/test/acos_r3fp64.f90
index 750eaa7..f91ca83 100644
--- a/test/acos_r3fp64.f90
+++ b/test/acos_r3fp64.f90
@@ -46,7 +46,7 @@ integer function acos_r3fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/acos_r4fp32.f90 b/test/acos_r4fp32.f90
index 24cfc7d..4c0da68 100644
--- a/test/acos_r4fp32.f90
+++ b/test/acos_r4fp32.f90
@@ -51,7 +51,7 @@ integer function acos_r4fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/acos_r4fp64.f90 b/test/acos_r4fp64.f90
index 0d7daef..3330ddc 100644
--- a/test/acos_r4fp64.f90
+++ b/test/acos_r4fp64.f90
@@ -51,7 +51,7 @@ integer function acos_r4fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/asin_r3fp32.f90 b/test/asin_r3fp32.f90
index 9e9cc88..18ad506 100644
--- a/test/asin_r3fp32.f90
+++ b/test/asin_r3fp32.f90
@@ -45,7 +45,7 @@ integer function asin_r3fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/asin_r3fp64.f90 b/test/asin_r3fp64.f90
index af99b68..5389322 100644
--- a/test/asin_r3fp64.f90
+++ b/test/asin_r3fp64.f90
@@ -45,7 +45,7 @@ integer function asin_r3fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/atan_r3fp32.f90 b/test/atan_r3fp32.f90
index 08705bf..7ada63f 100644
--- a/test/atan_r3fp32.f90
+++ b/test/atan_r3fp32.f90
@@ -46,7 +46,7 @@ integer function atan_r3fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/atan_r3fp64.f90 b/test/atan_r3fp64.f90
index 7140cec..1af50fb 100644
--- a/test/atan_r3fp64.f90
+++ b/test/atan_r3fp64.f90
@@ -46,7 +46,7 @@ integer function atan_r3fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/cos_r3fp32.f90 b/test/cos_r3fp32.f90
index 568915f..fca05da 100644
--- a/test/cos_r3fp32.f90
+++ b/test/cos_r3fp32.f90
@@ -46,7 +46,7 @@ integer function cos_r3fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/cos_r3fp64.f90 b/test/cos_r3fp64.f90
index 6a65b60..699038d 100644
--- a/test/cos_r3fp64.f90
+++ b/test/cos_r3fp64.f90
@@ -46,7 +46,7 @@ integer function cos_r3fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/linear_r3fp32.f90 b/test/linear_r3fp32.f90
index 8ba44aa..3e87eee 100644
--- a/test/linear_r3fp32.f90
+++ b/test/linear_r3fp32.f90
@@ -40,14 +40,14 @@ integer function linear_r3fp32() result(r)
   do k = 1,N
     do j = 1,N
       do i = 1,N
-        fexact(i,j,k) = (x(i,j,k,1)**3 - 1.0_real32)*(x(i,j,k,2)**3 - 1.0_real32)*(x(i,j,k,2)**3 - 1.0_real32)
+        fexact(i,j,k) = (x(i,j,k,1)**3 - 1.0_real32)*(x(i,j,k,2)**3 - 1.0_real32)*(x(i,j,k,3)**3 - 1.0_real32)
       end do
     end do
   end do
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= 10.0_real32*epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)*10.0_real32) then
     r = 0
   else
     r = 1
diff --git a/test/linear_r3fp64.f90 b/test/linear_r3fp64.f90
index aae28fe..98e4d9d 100644
--- a/test/linear_r3fp64.f90
+++ b/test/linear_r3fp64.f90
@@ -40,14 +40,14 @@ integer function linear_r3fp64() result(r)
   do k = 1,N
     do j = 1,N
       do i = 1,N
-        fexact(i,j,k) = (x(i,j,k,1)**3 - 1.0_real64)*(x(i,j,k,2)**3 - 1.0_real64)*(x(i,j,k,2)**3 - 1.0_real64)
+        fexact(i,j,k) = (x(i,j,k,1)**3 - 1.0_real64)*(x(i,j,k,2)**3 - 1.0_real64)*(x(i,j,k,3)**3 - 1.0_real64)
       end do
     end do
   end do
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= 10.0_real64*epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)*10.0_real64) then
     r = 0
   else
     r = 1
diff --git a/test/linear_r4fp32.f90 b/test/linear_r4fp32.f90
index f9f09c2..b906592 100644
--- a/test/linear_r4fp32.f90
+++ b/test/linear_r4fp32.f90
@@ -44,7 +44,7 @@ integer function linear_r4fp32() result(r)
     do k = 1,N
     do j = 1,N
       do i = 1,N
-        fexact(i,j,k,l) = (x(i,j,k,l,1)**3 - 1.0_real32)*(x(i,j,k,l,2)**3 - 1.0_real32)*(x(i,j,k,l,2)**3 - 1.0_real32)
+        fexact(i,j,k,l) = (x(i,j,k,l,1)**3 - 1.0_real32)*(x(i,j,k,l,2)**3 - 1.0_real32)*(x(i,j,k,l,3)**3 - 1.0_real32)
       end do
     end do
     end do
@@ -52,7 +52,7 @@ integer function linear_r4fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= 10.0_real32*epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)*10.0_real32) then
     r = 0
   else
     r = 1
diff --git a/test/linear_r4fp64.f90 b/test/linear_r4fp64.f90
index 1aa59b7..271e6f8 100644
--- a/test/linear_r4fp64.f90
+++ b/test/linear_r4fp64.f90
@@ -44,7 +44,7 @@ integer function linear_r4fp64() result(r)
     do k = 1,N
     do j = 1,N
       do i = 1,N
-        fexact(i,j,k,l) = (x(i,j,k,l,1)**3 - 1.0_real64)*(x(i,j,k,l,2)**3 - 1.0_real64)*(x(i,j,k,l,2)**3 - 1.0_real64)
+        fexact(i,j,k,l) = (x(i,j,k,l,1)**3 - 1.0_real64)*(x(i,j,k,l,2)**3 - 1.0_real64)*(x(i,j,k,l,3)**3 - 1.0_real64)
       end do
     end do
     end do
@@ -52,7 +52,7 @@ integer function linear_r4fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= 10.0_real64*epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)*10.0_real64) then
     r = 0
   else
     r = 1
diff --git a/test/log_r3fp32.f90 b/test/log_r3fp32.f90
index 58a1c93..1923504 100644
--- a/test/log_r3fp32.f90
+++ b/test/log_r3fp32.f90
@@ -38,14 +38,14 @@ integer function log_r3fp32() result(r)
   do k = 1,N
     do j = 1,N
       do i = 1,N
-        fexact(i,j,k) = log(x(i,j,k,1))*log(x(i,j,k,2))*log(x(i,j,k,2))
+        fexact(i,j,k) = log(x(i,j,k,1))*log(x(i,j,k,2))*log(x(i,j,k,3))
       end do
     end do
   end do
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/log_r3fp64.f90 b/test/log_r3fp64.f90
index cbeff2a..2fef1ab 100644
--- a/test/log_r3fp64.f90
+++ b/test/log_r3fp64.f90
@@ -38,14 +38,14 @@ integer function log_r3fp64() result(r)
   do k = 1,N
     do j = 1,N
       do i = 1,N
-        fexact(i,j,k) = log(x(i,j,k,1))*log(x(i,j,k,2))*log(x(i,j,k,2))
+        fexact(i,j,k) = log(x(i,j,k,1))*log(x(i,j,k,2))*log(x(i,j,k,3))
       end do
     end do
   end do
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/log_r4fp32.f90 b/test/log_r4fp32.f90
index 9f25582..8d83fad 100644
--- a/test/log_r4fp32.f90
+++ b/test/log_r4fp32.f90
@@ -42,7 +42,7 @@ integer function log_r4fp32() result(r)
     do k = 1,N
     do j = 1,N
       do i = 1,N
-        fexact(i,j,k,l) = log(x(i,j,k,l,1))*log(x(i,j,k,l,2))*log(x(i,j,k,l,2))
+        fexact(i,j,k,l) = log(x(i,j,k,l,1))*log(x(i,j,k,l,2))*log(x(i,j,k,l,3))
       end do
     end do
     end do
@@ -50,7 +50,7 @@ integer function log_r4fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/log_r4fp64.f90 b/test/log_r4fp64.f90
index 229adbe..5be8032 100644
--- a/test/log_r4fp64.f90
+++ b/test/log_r4fp64.f90
@@ -42,7 +42,7 @@ integer function log_r4fp64() result(r)
     do k = 1,N
     do j = 1,N
       do i = 1,N
-        fexact(i,j,k,l) = log(x(i,j,k,l,1))*log(x(i,j,k,l,2))*log(x(i,j,k,l,2))
+        fexact(i,j,k,l) = log(x(i,j,k,l,1))*log(x(i,j,k,l,2))*log(x(i,j,k,l,3))
       end do
     end do
     end do
@@ -50,7 +50,7 @@ integer function log_r4fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/sqrt_r3fp32.f90 b/test/sqrt_r3fp32.f90
index 21a5b6f..d21b106 100644
--- a/test/sqrt_r3fp32.f90
+++ b/test/sqrt_r3fp32.f90
@@ -45,7 +45,7 @@ integer function sqrt_r3fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/sqrt_r4fp32.f90 b/test/sqrt_r4fp32.f90
index 147d8b1..a5653c4 100644
--- a/test/sqrt_r4fp32.f90
+++ b/test/sqrt_r4fp32.f90
@@ -50,7 +50,7 @@ integer function sqrt_r4fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/sqrt_r4fp64.f90 b/test/sqrt_r4fp64.f90
index 14794e9..67d57ab 100644
--- a/test/sqrt_r4fp64.f90
+++ b/test/sqrt_r4fp64.f90
@@ -50,7 +50,7 @@ integer function sqrt_r4fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)*10.0_real64) then
     r = 0
   else
     r = 1
diff --git a/test/tan_r3fp32.f90 b/test/tan_r3fp32.f90
index dc5039c..34a43b9 100644
--- a/test/tan_r3fp32.f90
+++ b/test/tan_r3fp32.f90
@@ -46,7 +46,7 @@ integer function tan_r3fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/tan_r3fp64.f90 b/test/tan_r3fp64.f90
index fcb986f..bd96b2b 100644
--- a/test/tan_r3fp64.f90
+++ b/test/tan_r3fp64.f90
@@ -46,7 +46,7 @@ integer function tan_r3fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1
diff --git a/test/tan_r4fp32.f90 b/test/tan_r4fp32.f90
index 4dfcae5..3400500 100644
--- a/test/tan_r4fp32.f90
+++ b/test/tan_r4fp32.f90
@@ -51,7 +51,7 @@ integer function tan_r4fp32() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real32)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real32)) then
     r = 0
   else
     r = 1
diff --git a/test/tan_r4fp64.f90 b/test/tan_r4fp64.f90
index eb2be3c..a4df3ae 100644
--- a/test/tan_r4fp64.f90
+++ b/test/tan_r4fp64.f90
@@ -51,7 +51,7 @@ integer function tan_r4fp64() result(r)
 
   ! Evaluate the equation
   feval = f % evaluate(x)
-  if (maxval(abs(feval - fexact)) <= epsilon(1.0_real64)) then
+  if (maxval(abs(feval - fexact)) <= maxval(abs(fexact))*epsilon(1.0_real64)) then
     r = 0
   else
     r = 1