Browse Source

* first approach on better exception handling on hyperbolic harmonic functions, might need more fine-grained checking for good performance

florian 1 month ago
parent
commit
a9e735f257
1 changed files with 22 additions and 6 deletions
  1. 22 6
      rtl/objpas/math.pp

+ 22 - 6
rtl/objpas/math.pp

@@ -162,7 +162,7 @@ Const
        GreaterThanValue = types.GreaterThanValue;
        GreaterThanValue = types.GreaterThanValue;
 {$ENDIF FPC_DOTTEDUNITS}
 {$ENDIF FPC_DOTTEDUNITS}
 
 
-       
+
 {$push}
 {$push}
 {$R-}
 {$R-}
 {$Q-}
 {$Q-}
@@ -1290,10 +1290,11 @@ function cosh(x : Single) : Single;
   var
   var
      temp : ValReal;
      temp : ValReal;
   begin
   begin
+{$push}
+{$checkfpuexceptions on}
      if (x>8.94159862326326216608E+0001) or (x<-8.94159862326326216608E+0001) then
      if (x>8.94159862326326216608E+0001) or (x<-8.94159862326326216608E+0001) then
        exit(huge_single*huge_single);
        exit(huge_single*huge_single);
     temp:=exp(x);
     temp:=exp(x);
-{$push}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
      cosh:=0.5*(temp+1.0/temp);
      cosh:=0.5*(temp+1.0/temp);
 {$pop}
 {$pop}
@@ -1306,10 +1307,11 @@ function cosh(x : Double) : Double;
   var
   var
      temp : ValReal;
      temp : ValReal;
   begin
   begin
+{$push}
+{$checkfpuexceptions on}
      if (x>7.10475860073943942030E+0002) or (x<-7.10475860073943942030E+0002) then
      if (x>7.10475860073943942030E+0002) or (x<-7.10475860073943942030E+0002) then
        exit(huge_double*huge_double);
        exit(huge_double*huge_double);
      temp:=exp(x);
      temp:=exp(x);
-{$push}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
      cosh:=0.5*(temp+1.0/temp);
      cosh:=0.5*(temp+1.0/temp);
 {$pop}
 {$pop}
@@ -1333,6 +1335,8 @@ function sinh(x : Single) : Single;
   var
   var
      temp : ValReal;
      temp : ValReal;
   begin
   begin
+{$push}
+{$checkfpuexceptions on}
      if x>8.94159862326326216608E+0001 then
      if x>8.94159862326326216608E+0001 then
        exit(huge_single*huge_single);
        exit(huge_single*huge_single);
      if x<-8.94159862326326216608E+0001 then
      if x<-8.94159862326326216608E+0001 then
@@ -1341,7 +1345,6 @@ function sinh(x : Single) : Single;
      { gives better behavior around zero, and in particular ensures that sinh(-0.0)=-0.0 }
      { gives better behavior around zero, and in particular ensures that sinh(-0.0)=-0.0 }
      if temp=1 then
      if temp=1 then
        exit(x);
        exit(x);
-{$push}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
      sinh:=0.5*(temp-1.0/temp);
      sinh:=0.5*(temp-1.0/temp);
 {$pop}
 {$pop}
@@ -1354,6 +1357,8 @@ function sinh(x : Double) : Double;
   var
   var
      temp : ValReal;
      temp : ValReal;
   begin
   begin
+{$push}
+{$checkfpuexceptions on}
      if x>7.10475860073943942030E+0002 then
      if x>7.10475860073943942030E+0002 then
        exit(huge_double*huge_double);
        exit(huge_double*huge_double);
      if x<-7.10475860073943942030E+0002 then
      if x<-7.10475860073943942030E+0002 then
@@ -1361,7 +1366,6 @@ function sinh(x : Double) : Double;
      temp:=exp(x);
      temp:=exp(x);
      if temp=1 then
      if temp=1 then
        exit(x);
        exit(x);
-{$push}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
      sinh:=0.5*(temp-1.0/temp);
      sinh:=0.5*(temp-1.0/temp);
 {$pop}
 {$pop}
@@ -1392,6 +1396,7 @@ function tanh(x : Single) : Single;
       if tmp=1 then
       if tmp=1 then
         exit(x);
         exit(x);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
       result:=(tmp-1)/(1+tmp)
       result:=(tmp-1)/(1+tmp)
 {$pop}
 {$pop}
@@ -1401,6 +1406,7 @@ function tanh(x : Single) : Single;
       if tmp=1 then
       if tmp=1 then
         exit(x);
         exit(x);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
       result:=(1-tmp)/(1+tmp)
       result:=(1-tmp)/(1+tmp)
 {$pop}
 {$pop}
@@ -1419,6 +1425,7 @@ function tanh(x : Double) : Double;
       if tmp=1 then
       if tmp=1 then
         exit(x);
         exit(x);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
       result:=(tmp-1)/(1+tmp)
       result:=(tmp-1)/(1+tmp)
 {$pop}
 {$pop}
@@ -1428,6 +1435,7 @@ function tanh(x : Double) : Double;
       if tmp=1 then
       if tmp=1 then
         exit(x);
         exit(x);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
       result:=(1-tmp)/(1+tmp)
       result:=(1-tmp)/(1+tmp)
 {$pop}
 {$pop}
@@ -1466,6 +1474,7 @@ begin
   //SecH = 2 / (e^X + e^-X)
   //SecH = 2 / (e^X + e^-X)
   Ex:=Exp(X);
   Ex:=Exp(X);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
   SecH:=2/(Ex+1/Ex);
   SecH:=2/(Ex+1/Ex);
 {$pop}
 {$pop}
@@ -1480,6 +1489,7 @@ var
 begin
 begin
   Ex:=Exp(X);
   Ex:=Exp(X);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
   SecH:=2/(Ex+1/Ex);
   SecH:=2/(Ex+1/Ex);
 {$pop}
 {$pop}
@@ -1505,6 +1515,7 @@ begin
   //CscH = 2 / (e^X - e^-X)
   //CscH = 2 / (e^X - e^-X)
   Ex:=Exp(X);
   Ex:=Exp(X);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
   CscH:=2/(Ex-1/Ex);
   CscH:=2/(Ex-1/Ex);
 {$pop}
 {$pop}
@@ -1519,6 +1530,7 @@ var
 begin
 begin
   Ex:=Exp(X);
   Ex:=Exp(X);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
   CscH:=2/(Ex-1/Ex);
   CscH:=2/(Ex-1/Ex);
 {$pop}
 {$pop}
@@ -1546,6 +1558,7 @@ begin
     if e2=1 then
     if e2=1 then
       exit(1/x);
       exit(1/x);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
     result:=(1+e2)/(e2-1)
     result:=(1+e2)/(e2-1)
 {$pop}
 {$pop}
@@ -1555,6 +1568,7 @@ begin
     if e2=1 then
     if e2=1 then
       exit(1/x);
       exit(1/x);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
     result:=(1+e2)/(1-e2)
     result:=(1+e2)/(1-e2)
 {$pop}
 {$pop}
@@ -1573,6 +1587,7 @@ begin
     if e2=1 then
     if e2=1 then
       exit(1/x);
       exit(1/x);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
     result:=(1+e2)/(e2-1)
     result:=(1+e2)/(e2-1)
 {$pop}
 {$pop}
@@ -1582,6 +1597,7 @@ begin
     if e2=1 then
     if e2=1 then
       exit(1/x);
       exit(1/x);
 {$push}
 {$push}
+{$checkfpuexceptions on}
 {$safefpuexceptions on}
 {$safefpuexceptions on}
     result:=(1+e2)/(1-e2)
     result:=(1+e2)/(1-e2)
 {$pop}
 {$pop}
@@ -1959,7 +1975,7 @@ begin
       else
       else
         res:=0;
         res:=0;
       exit;
       exit;
-    end; 
+    end;
   res:=1;
   res:=1;
   while exponent<>0 do
   while exponent<>0 do
     begin
     begin