Changeset 12529

Timestamp:

06/25/12 12:16:37 (13 years ago)

Author:

Mathieu Morlighem

Message:

replaced all isnan by template xIsnan

Location:

issm/trunk-jpl/src/c

Files:

• modules/Bamgx/Bamgx.cpp (modified) (3 diffs)
• modules/Gradjx/Gradjx.cpp (modified) (1 diff)
• modules/InterpFromGridToMeshx/InterpFromGridToMeshx.cpp (modified) (3 diffs)
• modules/InterpFromMesh2dx/InterpFromMesh2dxt.cpp (modified) (1 diff)
• modules/InterpFromMeshToGridx/InterpFromMeshToGridx.cpp (modified) (1 diff)
• modules/InterpFromMeshToMesh3dx/InterpFromMeshToMesh3dx.cpp (modified) (1 diff)
• modules/IoModelToConstraintsx/IoModelToConstraintsx.cpp (modified) (2 diffs)
• modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp (modified) (12 diffs)
• modules/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp (modified) (1 diff)
• modules/ModelProcessorx/DiagnosticVert/CreateConstraintsDiagnosticVert.cpp (modified) (1 diff)
• modules/ModelProcessorx/Enthalpy/CreateConstraintsEnthalpy.cpp (modified) (1 diff)
• modules/ModelProcessorx/NodesPartitioning.cpp (modified) (1 diff)
• modules/ModelProcessorx/Thermal/CreateLoadsThermal.cpp (modified) (1 diff)
• objects/Bamg/Geometry.cpp (modified) (1 diff)
• objects/Constraints/SpcDynamic.cpp (modified) (1 diff)
• objects/Constraints/SpcStatic.cpp (modified) (1 diff)
• objects/Elements/Penta.cpp (modified) (17 diffs)
• objects/Elements/Tria.cpp (modified) (9 diffs)
• objects/Inputs/BoolInput.cpp (modified) (1 diff)
• objects/Inputs/DoubleInput.cpp (modified) (1 diff)
• objects/Inputs/IntInput.cpp (modified) (1 diff)
• objects/Inputs/PentaP1Input.cpp (modified) (1 diff)
• objects/Inputs/TriaP1Input.cpp (modified) (1 diff)
• objects/Loads/Friction.cpp (modified) (4 diffs)
• objects/Loads/Numericalflux.cpp (modified) (2 diffs)
• objects/Node.cpp (modified) (1 diff)
• objects/Numerics/ElementMatrix.cpp (modified) (1 diff)
• objects/Numerics/ElementVector.cpp (modified) (1 diff)
• shared/Numerics/BrentSearch.cpp (modified) (6 diffs)
• shared/Numerics/OptimalSearch.cpp (modified) (2 diffs)
• shared/Numerics/XZvectorsToCoordinateSystem.cpp (modified) (1 diff)
• shared/Numerics/isnan.cpp (modified) (1 diff)
• solutions/controlconvergence.cpp (modified) (1 diff)
• solutions/convergence.cpp (modified) (6 diffs)
• solutions/gradient_core.cpp (modified) (1 diff)

issm/trunk-jpl/src/c/modules/Bamgx/Bamgx.cpp

@@ -119 +119 @@
                         if (verbosity>1) _printLine_("   Merging Metric with hVertices...");
                         for (i=0;i<BTh.nbv;i++){
-                                if (!isnan(bamgopts->hVertices[i])){
+                                if (!xIsNan<IssmDouble>(bamgopts->hVertices[i])){
                                         BTh[i].m=Metric((float)bamgopts->hVertices[i]);
                                 }
@@ -129 +129 @@
                         if (verbosity>1) _printLine_("   Merging Metric with hminVertices...");
                         for (i=0;i<BTh.nbv;i++){
-                                if (!isnan(bamgopts->hminVertices[i])){
+                                if (!xIsNan<IssmDouble>(bamgopts->hminVertices[i])){
                                         Metric M=BTh.vertices[i].m;
                                         EigenMetric Vp(M/coef);
@@ -142 +142 @@
                         if (verbosity>1) _printLine_("   Merging Metric with hmaxVertices...");
                         for (i=0;i<BTh.nbv;i++){
-                                if (!isnan(bamgopts->hmaxVertices[i])){
+                                if (!xIsNan<IssmDouble>(bamgopts->hmaxVertices[i])){
                                         Metric M=BTh.vertices[i].m;
                                         EigenMetric Vp(M/coef);

issm/trunk-jpl/src/c/modules/Gradjx/Gradjx.cpp

@@ -55 +55 @@
         norm_inf=gradient->Norm(NORM_INF);
         if(norm_inf<=0)    _error2_("||∂J/∂α||∞ = 0    gradient norm is zero");
-        if(isnan(norm_inf))_error2_("||∂J/∂α||∞ = NaN  gradient norm is NaN");
+        if(xIsNan<IssmDouble>(norm_inf))_error2_("||∂J/∂α||∞ = NaN  gradient norm is NaN");
 
         /*Clean-up and assign output pointer*/

issm/trunk-jpl/src/c/modules/InterpFromGridToMeshx/InterpFromGridToMeshx.cpp

@@ -13 +13 @@
 #include "../../shared/shared.h"
 #include "../../include/include.h"
+#include "../../io/io.h"
 /*}}}*/
 
@@ -140 +141 @@
 
                 if(debug && my_thread==0)
-                 _pprintString_("\r      interpolation progress: "<<setw(6)<<setprecision(2)<<double(i-i0)/double(i1-i0)*100<<"%");
+                 _printString_("\r      interpolation progress: "<<setw(6)<<setprecision(2)<<double(i-i0)/double(i1-i0)*100<<"%");
                 x_grid=*(x_mesh+i);
                 y_grid=*(y_mesh+i);
@@ -179 +180 @@
                                         return NULL; /*WARNING: no error because it would blow up the multithreading!*/
                         }
-                        if(isnan(data_value)) data_value=default_value;
+                        if(xIsNan<IssmDouble>(data_value)) data_value=default_value;
                 }
                 else{

issm/trunk-jpl/src/c/modules/InterpFromMesh2dx/InterpFromMesh2dxt.cpp

@@ -113 +113 @@
                                                 data_value=data[i];
                                         }
-                                        if (isnan(data_value)){
+                                        if (xIsNan<IssmDouble>(data_value)){
                                                 if(num_default_values==1) data_value=default_values[0];
                                                 else data_value=default_values[j];

issm/trunk-jpl/src/c/modules/InterpFromMeshToGridx/InterpFromMeshToGridx.cpp

@@ -160 +160 @@
                                                 data_value=data_mesh[n];
                                         }
-                                        if (isnan(data_value)) data_value=default_value;
+                                        if (xIsNan<IssmDouble>(data_value)) data_value=default_value;
 
                                         /*insert value and go to the next point*/

issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh3dx/InterpFromMeshToMesh3dx.cpp

@@ -124 +124 @@
                                                 data_value=data[i];
                                         }
-                                        if (isnan(data_value)) data_value=default_value;
+                                        if (xIsNan<IssmDouble>(data_value)) data_value=default_value;
 
                                         /*insert value and go to the next point*/

issm/trunk-jpl/src/c/modules/IoModelToConstraintsx/IoModelToConstraintsx.cpp

@@ -55 +55 @@
                         if((iomodel->my_vertices[i])){
 
-                                if (!isnan(IssmDoublevector[i])){
+                                if (!xIsNan<IssmDouble>(IssmDoublevector[i])){
 
                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,IssmDoublevector[i],analysis_type));
@@ -87 +87 @@
                                 for(j=0;j<N;j++){
                                         values[j]=IssmDoublevector[i*N+j];
-                                        if(!isnan(values[j]))spcpresent=true; //NaN means no spc by default
+                                        if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
                                 }
 

issm/trunk-jpl/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp

@@ -120 +120 @@
                                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                 count++;
-                                                if (!isnan(spcvx[i])){
+                                                if (!xIsNan<IssmDouble>(spcvx[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
                                                 }
-                                                if (!isnan(spcvy[i])){
+                                                if (!xIsNan<IssmDouble>(spcvy[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
@@ -135 +135 @@
                                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                 count++;
-                                                if (!isnan(spcvx[i])){
+                                                if (!xIsNan<IssmDouble>(spcvx[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
                                                 }
-                                                if (!isnan(spcvy[i])){
+                                                if (!xIsNan<IssmDouble>(spcvy[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
@@ -157 +157 @@
                                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                 count++;
-                                                if (!isnan(spcvx[i])){
+                                                if (!xIsNan<IssmDouble>(spcvx[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
                                                 }
-                                                if (!isnan(spcvy[i])){
+                                                if (!xIsNan<IssmDouble>(spcvy[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
@@ -172 +172 @@
                                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                 count++;
-                                                if (!isnan(spcvx[i])){
+                                                if (!xIsNan<IssmDouble>(spcvx[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
                                                 }
-                                                if (!isnan(spcvy[i])){
+                                                if (!xIsNan<IssmDouble>(spcvy[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
                                                 }
-                                                if (!isnan(spcvz[i])){
+                                                if (!xIsNan<IssmDouble>(spcvz[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,spcvz[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
@@ -197 +197 @@
                                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                 count++;
-                                                if (!isnan(spcvx[i])){
+                                                if (!xIsNan<IssmDouble>(spcvx[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
                                                 }
-                                                if (!isnan(spcvy[i])){
+                                                if (!xIsNan<IssmDouble>(spcvy[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
@@ -212 +212 @@
                                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,0,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                 count++;
-                                                if (!isnan(spcvx[i])){
+                                                if (!xIsNan<IssmDouble>(spcvx[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
                                                 }
-                                                if (!isnan(spcvy[i])){
+                                                if (!xIsNan<IssmDouble>(spcvy[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
                                                 }
-                                                if (!isnan(spcvz[i])){
+                                                if (!xIsNan<IssmDouble>(spcvz[i])){
                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,spcvz[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                                         count++;
@@ -229 +229 @@
                         /*Now add the regular spcs*/
                         else{
-                                if (Mx==numberofvertices && !isnan(spcvx[i])){
+                                if (Mx==numberofvertices && !xIsNan<IssmDouble>(spcvx[i])){
                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                         count++;
@@ -239 +239 @@
                                         for(j=0;j<Nx;j++){
                                                 values[j]=spcvx[i*Nx+j]/yts;
-                                                if(!isnan(values[j]))spcpresent=true; //NaN means no spc by default
+                                                if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
                                         }
 
@@ -253 +253 @@
                                 }
 
-                                if (My==numberofvertices && !isnan(spcvy[i])){
+                                if (My==numberofvertices && !xIsNan<IssmDouble>(spcvy[i])){
                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vy.
                                         count++;
@@ -263 +263 @@
                                         for(j=0;j<Ny;j++){
                                                 values[j]=spcvy[i*Ny+j]/yts;
-                                                if(!isnan(values[j]))spcpresent=true; //NaN means no spc by default
+                                                if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
                                         }
                                         if(spcpresent){
@@ -277 +277 @@
 
                                 if ((int)vertices_type[i]==StokesApproximationEnum ||  ((int)vertices_type[i]==NoneApproximationEnum)){
-                                        if (Mz==numberofvertices && !isnan(spcvz[i])){
+                                        if (Mz==numberofvertices && !xIsNan<IssmDouble>(spcvz[i])){
                                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,spcvz[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy
                                                 count++;
@@ -287 +287 @@
                                                 for(j=0;j<Nz;j++){
                                                         values[j]=spcvz[i*Nz+j]/yts;
-                                                        if(!isnan(values[j]))spcpresent=true; //NaN means no spc by default
+                                                        if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
                                                 }
                                                 if(spcpresent){

issm/trunk-jpl/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp

@@ -60 +60 @@
                         }
                         else{
-                                if (!isnan(iomodel->Data(DiagnosticSpcvxEnum)[i])){
+                                if (!xIsNan<IssmDouble>(iomodel->Data(DiagnosticSpcvxEnum)[i])){
                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->Data(DiagnosticSpcvxEnum)[i]/yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
                                         count++;
                                 }
 
-                                if (!isnan(iomodel->Data(DiagnosticSpcvyEnum)[i])){
+                                if (!xIsNan<IssmDouble>(iomodel->Data(DiagnosticSpcvyEnum)[i])){
                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->Data(DiagnosticSpcvyEnum)[i]/yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy
                                         count++;

issm/trunk-jpl/src/c/modules/ModelProcessorx/DiagnosticVert/CreateConstraintsDiagnosticVert.cpp

@@ -56 +56 @@
                                 count++;
                         }
-                        else if (!isnan(iomodel->Data(DiagnosticSpcvzEnum)[i])){
+                        else if (!xIsNan<IssmDouble>(iomodel->Data(DiagnosticSpcvzEnum)[i])){
                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,
                                                                 iomodel->Data(DiagnosticSpcvzEnum)[i]/yts,DiagnosticVertAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.

issm/trunk-jpl/src/c/modules/ModelProcessorx/Enthalpy/CreateConstraintsEnthalpy.cpp

@@ -55 +55 @@
                 if((iomodel->my_vertices[i])){
 
-                        if (!isnan(spctemperature[i])){
+                        if (!xIsNan<IssmDouble>(spctemperature[i])){
 
                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,heatcapacity*(spctemperature[i]-referencetemperature),EnthalpyAnalysisEnum));

issm/trunk-jpl/src/c/modules/ModelProcessorx/NodesPartitioning.cpp

@@ -125 +125 @@
                  * we must clone the nodes on this partition so that the loads (Numericalflux)
                  * will have access to their properties (dofs,...)*/
-                if(my_elements[(int)e1] && !isnan(e2) && !my_elements[(int)e2]){ 
+                if(my_elements[(int)e1] && !xIsNan<IssmDouble>(e2) && !my_elements[(int)e2]){ 
 
                         /*1: Get vertices ids*/

issm/trunk-jpl/src/c/modules/ModelProcessorx/Thermal/CreateLoadsThermal.cpp

@@ -43 +43 @@
                 /*keep only this partition's nodes:*/
                 if((iomodel->my_vertices[i]==1)){
-                        if (isnan(iomodel->Data(ThermalSpctemperatureEnum)[i])){ //No penalty applied on spc nodes!
+                        if (xIsNan<IssmDouble>(iomodel->Data(ThermalSpctemperatureEnum)[i])){ //No penalty applied on spc nodes!
                                 loads->AddObject(new Pengrid(iomodel->loadcounter+i+1,i,iomodel,ThermalAnalysisEnum));
                         }

issm/trunk-jpl/src/c/objects/Bamg/Geometry.cpp

@@ -183 +183 @@
                         if(verbose>5) _printLine_("      processing hVertices");
                         for (i=0;i< nbv;i++){
-                                if (!isnan(bamgopts->hVertices[i])){
+                                if (!xIsNan<IssmDouble>(bamgopts->hVertices[i])){
                                         vertices[i].m=Metric((double)bamgopts->hVertices[i]);
                                 }

issm/trunk-jpl/src/c/objects/Constraints/SpcDynamic.cpp

@@ -124 +124 @@
 IssmDouble SpcDynamic::GetValue(){
         _assert_(this->isset);
-        _assert_(!isnan(value));
+        _assert_(!xIsNan<IssmDouble>(value));
         return value;
 }

issm/trunk-jpl/src/c/objects/Constraints/SpcStatic.cpp

@@ -124 +124 @@
 /*FUNCTION SpcStatic::GetValue {{{*/
 IssmDouble SpcStatic::GetValue(){
-        _assert_(!isnan(value));
+        _assert_(!xIsNan<IssmDouble>(value));
         return value;
 }

issm/trunk-jpl/src/c/objects/Elements/Penta.cpp

@@ -66 +66 @@
 
         /*Build neighbors list*/
-        if (isnan(iomodel->Data(MeshUpperelementsEnum)[index])) penta_elements_ids[1]=this->id; //upper penta is the same penta
+        if (xIsNan<IssmDouble>(iomodel->Data(MeshUpperelementsEnum)[index])) penta_elements_ids[1]=this->id; //upper penta is the same penta
         else                                    penta_elements_ids[1]=(int)(iomodel->Data(MeshUpperelementsEnum)[index]);
-        if (isnan(iomodel->Data(MeshLowerelementsEnum)[index])) penta_elements_ids[0]=this->id; //lower penta is the same penta
+        if (xIsNan<IssmDouble>(iomodel->Data(MeshLowerelementsEnum)[index])) penta_elements_ids[0]=this->id; //lower penta is the same penta
         else                                    penta_elements_ids[0]=(int)(iomodel->Data(MeshLowerelementsEnum)[index]);
         this->InitHookNeighbors(penta_elements_ids);
@@ -1742 +1742 @@
         for(i=0;i<numdof2d;i++){
                 newthickness[i]=solution[doflist[i]];
-                if(isnan(newthickness[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(newthickness[i])) _error2_("NaN found in solution vector");
                 /*Constrain thickness to be at least 1m*/
                 if(newthickness[i]<minthickness) newthickness[i]=minthickness;
@@ -1812 +1812 @@
         for(int i=0;i<numdof;i++){
                 values[i]=solution[doflist[i]];
-                if(isnan(values[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(values[i])) _error2_("NaN found in solution vector");
         }
 
@@ -1842 +1842 @@
                 values[i]         =solution[doflist[i]];
                 values[i+numdof2d]=values[i];
-                if(isnan(values[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(values[i])) _error2_("NaN found in solution vector");
         }
 
@@ -4328 +4328 @@
 
                 /*Check solution*/
-                if(isnan(values[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(values[i])) _error2_("NaN found in solution vector");
                 //if(values[i]<0)      _printLine_("temperature < 0°K found in solution vector");
                 //if(values[i]>275)    _printLine_("temperature > 275°K found in solution vector (Paterson's rheology associated is negative)");
@@ -4398 +4398 @@
 
                 /*Check solution*/
-                if(isnan(values[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(values[i])) _error2_("NaN found in solution vector");
         }
 
@@ -4943 +4943 @@
                 /*Add gradje_g_gaussian vector to gradje_g: */
                 for(i=0;i<NUMVERTICES;i++){
-                        _assert_(!isnan(grade_g[i]));
+                        _assert_(!xIsNan<IssmDouble>(grade_g[i]));
                         grade_g[i]+=grade_g_gaussian[i];
                 }
@@ -5163 +5163 @@
 
                 /*Check solution*/
-                if(isnan(lambdax[i])) _error2_("NaN found in solution vector");
-                if(isnan(lambday[i])) _error2_("NaN found in solution vector");
-                if(isnan(lambdaz[i])) _error2_("NaN found in solution vector");
-                if(isnan(lambdap[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambdax[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambday[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambdaz[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambdap[i])) _error2_("NaN found in solution vector");
         }
 
@@ -5202 +5202 @@
 
                 /*Check solution*/
-                if(isnan(lambdax[i]))       _error2_("NaN found in solution vector");
-                if(isnan(lambday[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambdax[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambday[i]))       _error2_("NaN found in solution vector");
         }
 
@@ -8128 +8128 @@
 
                 /*Check solution*/
-                if(isnan(vx[i])) _error2_("NaN found in solution vector");
-                if(isnan(vy[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i])) _error2_("NaN found in solution vector");
         }
 
@@ -8226 +8226 @@
 
                 /*Check solution*/
-                if(isnan(vx[i])) _error2_("NaN found in solution vector");
-                if(isnan(vy[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i])) _error2_("NaN found in solution vector");
         }
 
@@ -8314 +8314 @@
 
                 /*Check solution*/
-                if(isnan(vx[i]))       _error2_("NaN found in solution vector");
-                if(isnan(vy[i]))       _error2_("NaN found in solution vector");
-                if(isnan(vzstokes[i])) _error2_("NaN found in solution vector");
-                if(isnan(pressure[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vzstokes[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(pressure[i])) _error2_("NaN found in solution vector");
         }
 
@@ -8393 +8393 @@
 
                 /*Check solution*/
-                if(isnan(vx[i])) _error2_("NaN found in solution vector");
-                if(isnan(vy[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i])) _error2_("NaN found in solution vector");
         }
 
@@ -8482 +8482 @@
 
                 /*Check solution*/
-                if(isnan(vx[i]))       _error2_("NaN found in solution vector");
-                if(isnan(vy[i]))       _error2_("NaN found in solution vector");
-                if(isnan(vzstokes[i])) _error2_("NaN found in solution vector");
-                if(isnan(pressure[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vzstokes[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(pressure[i])) _error2_("NaN found in solution vector");
         }
 
@@ -8558 +8558 @@
 
                 /*Check solution*/
-                if(isnan(vx[i])) _error2_("NaN found in solution vector");
-                if(isnan(vy[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i])) _error2_("NaN found in solution vector");
         }
 
@@ -8627 +8627 @@
 
                 /*Check solution*/
-                if(isnan(vz[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vz[i])) _error2_("NaN found in solution vector");
         }
 
@@ -8725 +8725 @@
 
                 /*Check solution*/
-                if(isnan(vx[i]))       _error2_("NaN found in solution vector");
-                if(isnan(vy[i]))       _error2_("NaN found in solution vector");
-                if(isnan(vz[i]))       _error2_("NaN found in solution vector");
-                if(isnan(pressure[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vz[i]))       _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(pressure[i])) _error2_("NaN found in solution vector");
         }
 

issm/trunk-jpl/src/c/objects/Elements/Tria.cpp

@@ -1523 +1523 @@
         for(int i=0;i<numdof;i++){
                 values[i]=solution[doflist[i]];
-                if(isnan(values[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(values[i])) _error2_("NaN found in solution vector");
         }
 
@@ -1556 +1556 @@
         for(i=0;i<numdof;i++){
                 newthickness[i]=solution[doflist[i]];
-                if(isnan(newthickness[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(newthickness[i])) _error2_("NaN found in solution vector");
                 /*Constrain thickness to be at least 1m*/
                 if(newthickness[i]<minthickness) newthickness[i]=minthickness;
@@ -3325 +3325 @@
 
                 /*Check solution*/
-                if(isnan(vx[i])) _error2_("NaN found in solution vector");
-                if(isnan(vy[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i])) _error2_("NaN found in solution vector");
         }
 
@@ -3385 +3385 @@
 
                 /*Check solution*/
-                if(isnan(vx[i])) _error2_("NaN found in solution vector");
-                if(isnan(vy[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vx[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(vy[i])) _error2_("NaN found in solution vector");
         }
 
@@ -3731 +3731 @@
                 /*Add gradje_g_gaussian vector to gradje_g: */
                 for(i=0;i<NUMVERTICES;i++){
-                        _assert_(!isnan(grade_g[i]));
+                        _assert_(!xIsNan<IssmDouble>(grade_g[i]));
                         grade_g[i]+=grade_g_gaussian[i];
                 }
@@ -3792 +3792 @@
                 for (i=0;i<NUMVERTICES;i++){
                         grade_g[i]+=-weight*Jdet*gauss->weight*(dbasis[0][i]*dk[0]+dbasis[1][i]*dk[1]);
-                        _assert_(!isnan(grade_g[i]));
+                        _assert_(!xIsNan<IssmDouble>(grade_g[i]));
                 }
         }
@@ -4947 +4947 @@
 
                 /*Check solution*/
-                if(isnan(lambdax[i])) _error2_("NaN found in solution vector");
-                if(isnan(lambday[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambdax[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambday[i])) _error2_("NaN found in solution vector");
         }
 
@@ -4978 +4978 @@
         for(i=0;i<numdof;i++){
                 lambda[i]=values[i];
-                if(isnan(lambda[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(lambda[i])) _error2_("NaN found in solution vector");
         }
 
@@ -5314 +5314 @@
         for(i=0;i<numdof;i++){
                 values[i]=solution[doflist[i]];
-                if(isnan(values[i])) _error2_("NaN found in solution vector");
+                if(xIsNan<IssmDouble>(values[i])) _error2_("NaN found in solution vector");
                 if (values[i]<pow((IssmDouble)10,(IssmDouble)-10))values[i]=pow((IssmDouble)10,(IssmDouble)-10); //correcting the water column to positive values
  

issm/trunk-jpl/src/c/objects/Inputs/BoolInput.cpp

@@ -169 +169 @@
 void BoolInput::Constrain(IssmPDouble cm_min, IssmPDouble cm_max){
 
-        if(!isnan(cm_min)) if (this->value<cm_min)this->value=cm_min;
-        if(!isnan(cm_max)) if (this->value>cm_max)this->value=cm_max;
+        if(!xIsNan<IssmDouble>(cm_min)) if (this->value<cm_min)this->value=cm_min;
+        if(!xIsNan<IssmDouble>(cm_max)) if (this->value>cm_max)this->value=cm_max;
 
 }

issm/trunk-jpl/src/c/objects/Inputs/DoubleInput.cpp

@@ -225 +225 @@
 void DoubleInput::Constrain(IssmPDouble cm_min, IssmPDouble cm_max){
 
-        if(!isnan(cm_min)) if (this->value<cm_min)this->value=cm_min;
-        if(!isnan(cm_max)) if (this->value>cm_max)this->value=cm_max;
+        if(!xIsNan<IssmDouble>(cm_min)) if (this->value<cm_min)this->value=cm_min;
+        if(!xIsNan<IssmDouble>(cm_max)) if (this->value>cm_max)this->value=cm_max;
 
 }

issm/trunk-jpl/src/c/objects/Inputs/IntInput.cpp

@@ -174 +174 @@
 void IntInput::Constrain(IssmPDouble cm_min, IssmPDouble cm_max){
 
-        if(!isnan(cm_min)) if (this->value<cm_min)this->value=(int)cm_min;
-        if(!isnan(cm_max)) if (this->value>cm_max)this->value=(int)cm_max;
+        if(!xIsNan<IssmDouble>(cm_min)) if (this->value<cm_min)this->value=(int)cm_min;
+        if(!xIsNan<IssmDouble>(cm_max)) if (this->value>cm_max)this->value=(int)cm_max;
 
 }

issm/trunk-jpl/src/c/objects/Inputs/PentaP1Input.cpp

@@ -471 +471 @@
         const int numnodes=6;
                 
-        if(!isnan(cm_min)) for(i=0;i<numnodes;i++)if (this->values[i]<cm_min)this->values[i]=cm_min;
-        if(!isnan(cm_max)) for(i=0;i<numnodes;i++)if (this->values[i]>cm_max)this->values[i]=cm_max;
+        if(!xIsNan<IssmDouble>(cm_min)) for(i=0;i<numnodes;i++)if (this->values[i]<cm_min)this->values[i]=cm_min;
+        if(!xIsNan<IssmDouble>(cm_max)) for(i=0;i<numnodes;i++)if (this->values[i]>cm_max)this->values[i]=cm_max;
 
 }

issm/trunk-jpl/src/c/objects/Inputs/TriaP1Input.cpp

@@ -334 +334 @@
         const int numnodes=3;
                 
-        if(!isnan(cm_min)) for(i=0;i<numnodes;i++)if (this->values[i]<cm_min)this->values[i]=cm_min;
-        if(!isnan(cm_max)) for(i=0;i<numnodes;i++)if (this->values[i]>cm_max)this->values[i]=cm_max;
+        if(!xIsNan<IssmDouble>(cm_min)) for(i=0;i<numnodes;i++)if (this->values[i]<cm_min)this->values[i]=cm_min;
+        if(!xIsNan<IssmDouble>(cm_max)) for(i=0;i<numnodes;i++)if (this->values[i]>cm_max)this->values[i]=cm_max;
 
 }

issm/trunk-jpl/src/c/objects/Loads/Friction.cpp

@@ -113 +113 @@
 
         alpha2=pow(drag_coefficient,2)*pow(Neff,r)*pow(vmag,(s-1));
-        _assert_(!isnan(alpha2));
+        _assert_(!xIsNan<IssmDouble>(alpha2));
 
         /*Assign output pointers:*/
@@ -177 +177 @@
 
         alpha2=pow(drag_coefficient,2)*pow(Neff,r)*pow(vmag,(s-1));
-        _assert_(!isnan(alpha2));
+        _assert_(!xIsNan<IssmDouble>(alpha2));
 
         /*Assign output pointers:*/
@@ -243 +243 @@
         if(vmag==0 && (s-1)<0) _error2_("velocity is 0 and (s-1)=" << (s-1) << "<0, alpha_complement is Inf");
 
-        alpha_complement=pow(Neff,r)*pow(vmag,(s-1));            _assert_(!isnan(alpha_complement));
+        alpha_complement=pow(Neff,r)*pow(vmag,(s-1));            _assert_(!xIsNan<IssmDouble>(alpha_complement));
 
         /*Assign output pointers:*/
@@ -309 +309 @@
         if(vmag==0 && (s-1)<0) _error2_("velocity is 0 and (s-1)=" << (s-1) << "<0, alpha_complement is Inf");
 
-        alpha_complement=pow(Neff,r)*pow(vmag,(s-1));            _assert_(!isnan(alpha_complement));
+        alpha_complement=pow(Neff,r)*pow(vmag,(s-1));            _assert_(!xIsNan<IssmDouble>(alpha_complement));
 
         /*Assign output pointers:*/

issm/trunk-jpl/src/c/objects/Loads/Numericalflux.cpp

@@ -60 +60 @@
 
         /*First, see wether this is an internal or boundary edge (if e2=NaN)*/
-        if (isnan((IssmDouble)iomodel->Data(MeshEdgesEnum)[4*i+3])){ //edges are [node1 node2 elem1 elem2]
+        if (xIsNan<IssmDouble>((IssmDouble)iomodel->Data(MeshEdgesEnum)[4*i+3])){ //edges are [node1 node2 elem1 elem2]
                 /* Boundary edge, only one element */
                 e1=(int)iomodel->Data(MeshEdgesEnum)[4*i+2];
@@ -751 +751 @@
                 vyaverage_input->GetInputValue(&vy,gauss);
                 spcthickness_input->GetInputValue(&thickness,gauss);
-                if(isnan(thickness)) _error2_("Cannot weakly apply constraint because NaN was provided");
+                if(xIsNan<IssmDouble>(thickness)) _error2_("Cannot weakly apply constraint because NaN was provided");
 
                 UdotN=vx*normal[0]+vy*normal[1];

issm/trunk-jpl/src/c/objects/Node.cpp

@@ -524 +524 @@
                         if(this->indexing.s_set[i]){
                                 values[count]=this->indexing.svalues[i];
-                                _assert_(!isnan(values[count]));
+                                _assert_(!xIsNan<IssmDouble>(values[count]));
                                 count++;
                         }

issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.cpp

@@ -343 +343 @@
         for (int i=0;i<this->nrows;i++){
                 for(int j=0;j<this->ncols;j++){
-                        if (isnan(this->values[i*this->ncols+j])) _error2_("NaN found in Element Matrix");
+                        if (xIsNan<IssmDouble>(this->values[i*this->ncols+j])) _error2_("NaN found in Element Matrix");
                         if (fabs(this->values[i*this->ncols+j])>1.e+50) _error2_("Element Matrix values exceeds 1.e+50");
                 }

issm/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp

@@ -210 +210 @@
 #ifdef _ISSM_DEBUG_ 
         for (int i=0;i<this->nrows;i++){
-                if (isnan(this->values[i])) _error2_("NaN found in Element Vector");
+                if (xIsNan<IssmDouble>(this->values[i])) _error2_("NaN found in Element Vector");
                 if (fabs( this->values[i])>1.e+50) _error2_("Element Vector values exceeds 1.e+50");
         }

issm/trunk-jpl/src/c/shared/Numerics/BrentSearch.cpp

@@ -44 +44 @@
         iter=0;
         fxmin = (*f)(xmin,optargs);
-        if (isnan(fxmin)) _error2_("Function evaluation returned NaN");
+        if (xIsNan<IssmDouble>(fxmin)) _error2_("Function evaluation returned NaN");
         cout<<setprecision(5);
         if(VerboseControl()) _pprintLine_("");
@@ -51 +51 @@
         if(VerboseControl()) _pprintLine_("           N/A    "<<setw(12)<<xmin<<"  "<<setw(12)<<fxmin<<"           N/A         boundary");
         fxmax = (*f)(xmax,optargs);
-        if (isnan(fxmax)) _error2_("Function evaluation returned NaN");
+        if (xIsNan<IssmDouble>(fxmax)) _error2_("Function evaluation returned NaN");
         if(VerboseControl()) _pprintLine_("           N/A    "<<setw(12)<<xmax<<"  "<<setw(12)<<fxmax<<"           N/A         boundary");
 
         /*test if jump option activated and xmin==0*/
-        if (!isnan(cm_jump) && (xmin==0) && (fxmax/fxmin)<cm_jump){
+        if (!xIsNan<IssmDouble>(cm_jump) && (xmin==0) && (fxmax/fxmin)<cm_jump){
                 *psearch_scalar=xmax;
                 *pJ=fxmax;
@@ -75 +75 @@
         /*2: call the function to be evaluated*/
         fxbest = (*f)(x,optargs);
-        if(isnan(fxbest)) _error2_("Function evaluation returned NaN");
+        if(xIsNan<IssmDouble>(fxbest)) _error2_("Function evaluation returned NaN");
         iter=iter+1;
 
@@ -90 +90 @@
         if(VerboseControl())
          _pprintLine_("         "<<setw(5)<<iter<<"    "<<setw(12)<<xbest<<"  "<<setw(12)<<fxbest<<"  "<<setw(12)<<pow(pow(xbest-xm,2),0.5)<<"         initial");
-        if (!isnan(cm_jump) && (xmin==0) && ((fxbest/fxmin)<cm_jump)){
+        if (!xIsNan<IssmDouble>(cm_jump) && (xmin==0) && ((fxbest/fxmin)<cm_jump)){
                 if(VerboseControl()) _pprintLine_("      optimization terminated: current x satisfies criteria 'cm_jump'=" << cm_jump);
                 loop=false;
@@ -159 +159 @@
                 //evaluate function on x
                 fx = (*f)(x,optargs);
-                if(isnan(fx)) _error2_("Function evaluation returned NaN");
+                if(xIsNan<IssmDouble>(fx)) _error2_("Function evaluation returned NaN");
                 iter=iter+1;
 
@@ -197 +197 @@
                         loop=false;
                 }
-                else if (!isnan(cm_jump) && (xmin==0) && ((fxbest/fxmin)<cm_jump)){
+                else if (!xIsNan<IssmDouble>(cm_jump) && (xmin==0) && ((fxbest/fxmin)<cm_jump)){
                         if(VerboseControl()) _pprintLine_("      optimization terminated: current x satisfies criteria 'cm_jump'=" << cm_jump);
                         loop=false;

issm/trunk-jpl/src/c/shared/Numerics/OptimalSearch.cpp

@@ -41 +41 @@
         //get the value of the function at the first boundary
         fx1= (*f)(x1,optargs);
-        if (isnan(fx1)) _error2_("Function evaluation returned NaN");
+        if (xIsNan<IssmDouble>(fx1)) _error2_("Function evaluation returned NaN");
         cout<<setprecision(5);
         if(VerboseControl()) _pprintLine_("");
@@ -57 +57 @@
                 iter++;
                 fx2 = (*f)(x2,optargs);
-                if (isnan(fx2)) _error2_("Function evaluation returned NaN");
+                if (xIsNan<IssmDouble>(fx2)) _error2_("Function evaluation returned NaN");
                 if(VerboseControl())
                  _pprintLine_("         "<<setw(5)<<iter<<"    "<<setw(12)<<x2<<"  "<<setw(12)<<fx2<<"  "<<(fabs(x2-x1)>fabs(fx2-fx1)?fabs(fx2-fx1):fabs(x2-x1)));

issm/trunk-jpl/src/c/shared/Numerics/XZvectorsToCoordinateSystem.cpp

@@ -12 +12 @@
 
         for(i=0;i<6;i++){
-                if(isnan(xzvectors[i])){
+                if(xIsNan<IssmDouble>(xzvectors[i])){
                         /*At least one NaN found: default to Id*/
                         T[0*3+0] = 1.0; T[0*3+1] = 0.0; T[0*3+2] = 0.0;

issm/trunk-jpl/src/c/shared/Numerics/isnan.cpp

@@ -2 +2 @@
 #ifdef _INTEL_WIN_
 
-int isnan(IssmPDouble x){
+int xIsNan<IssmDouble>(IssmPDouble x){
         return (x!=x)?1:0;
 }

issm/trunk-jpl/src/c/solutions/controlconvergence.cpp

@@ -23 +23 @@
 
         /*Has convergence been reached?*/
-        if (!isnan(tol_cm) && J<tol_cm){
+        if (!xIsNan<IssmDouble>(tol_cm) && J<tol_cm){
                 converged=true;
                 if(VerboseConvergence()) _pprintString_("      Convergence criterion reached: J = " << J << " < " << tol_cm);

issm/trunk-jpl/src/c/solutions/convergence.cpp

@@ -70 +70 @@
         nF=pf->Norm(NORM_TWO);
         res=nKUoldF/nF;
-        if (isnan(res)){
+        if (xIsNan<IssmDouble>(res)){
                 _pprintLine_("norm nf = " << nF << "f and norm kuold = " << nKUoldF << "f");
                 _error2_("mechanical equilibrium convergence criterion is NaN!");
@@ -90 +90 @@
 
         /*Relative criterion (optional)*/
-        if (!isnan(eps_rel) || (VerboseConvergence())){
+        if (!xIsNan<IssmDouble>(eps_rel) || (VerboseConvergence())){
 
                 //compute norm(du)/norm(u)
@@ -96 +96 @@
                 ndu=duf->Norm(NORM_TWO); nu=old_uf->Norm(NORM_TWO);
 
-                if (isnan(ndu) || isnan(nu)) _error2_("convergence criterion is NaN!");
+                if (xIsNan<IssmDouble>(ndu) || xIsNan<IssmDouble>(nu)) _error2_("convergence criterion is NaN!");
 
                 //clean up
@@ -102 +102 @@
 
                 //print
-                if (!isnan(eps_rel)){
+                if (!xIsNan<IssmDouble>(eps_rel)){
                         if((ndu/nu)<eps_rel){
                                 if(VerboseConvergence()) _pprintLine_("   Convergence criterion: norm(du)/norm(u)" << "50s" << ndu/nu*100 << " < " << eps_rel*100 << " %");
@@ -116 +116 @@
 
         /*Absolute criterion (Optional) = max(du)*/
-        if (!isnan(eps_abs) || (VerboseConvergence())){
+        if (!xIsNan<IssmDouble>(eps_abs) || (VerboseConvergence())){
 
                 //compute max(du)
                 duf=old_uf->Duplicate(); old_uf->Copy(duf); duf->AYPX(uf,-1.0);
                 ndu=duf->Norm(NORM_TWO); nduinf=duf->Norm(NORM_INF);
-                if (isnan(ndu) || isnan(nu)) _error2_("convergence criterion is NaN!");
+                if (xIsNan<IssmDouble>(ndu) || xIsNan<IssmDouble>(nu)) _error2_("convergence criterion is NaN!");
 
                 //clean up
@@ -127 +127 @@
 
                 //print
-                if (!isnan(eps_abs)){
+                if (!xIsNan<IssmDouble>(eps_abs)){
                         if ((nduinf*yts)<eps_abs){
                                 if(VerboseConvergence()) _pprintLine_("   Convergence criterion: max(du)" << "50s" << nduinf*yts << " < " << eps_abs << " m/yr");

issm/trunk-jpl/src/c/solutions/gradient_core.cpp

@@ -38 +38 @@
         norm_inf=new_gradient->Norm(NORM_INF);
         if(norm_inf<=0)    _error2_("||∂J/∂α||∞ = 0    gradient norm is zero");
-        if(isnan(norm_inf))_error2_("||∂J/∂α||∞ = NaN  gradient norm is NaN");
+        if(xIsNan<IssmDouble>(norm_inf))_error2_("||∂J/∂α||∞ = NaN  gradient norm is NaN");
 
         /*plug back into inputs: */