Changeset 20669


Ignore:
Timestamp:
05/30/16 12:40:38 (9 years ago)
Author:
Mathieu Morlighem
Message:

NEW: added check for isinf in solution vector

Location:
issm/trunk-jpl/src/c
Files:
19 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/analyses/AdjointHorizAnalysis.cpp

    r20608 r20669  
    23432343        /*fill in all arrays: */
    23442344        for(i=0;i<vnumnodes;i++){
    2345                 lambdax[i] = values[i*dim+0]; if(xIsNan<IssmDouble>(lambdax[i])) _error_("NaN found in solution vector");
    2346                 lambday[i] = values[i*dim+1]; if(xIsNan<IssmDouble>(lambday[i])) _error_("NaN found in solution vector");
     2345                lambdax[i] = values[i*dim+0];
     2346                if(xIsNan<IssmDouble>(lambdax[i])) _error_("NaN found in solution vector");
     2347                if(xIsInf<IssmDouble>(lambdax[i])) _error_("Inf found in solution vector");
     2348                lambday[i] = values[i*dim+1];
     2349                if(xIsNan<IssmDouble>(lambday[i])) _error_("NaN found in solution vector");
     2350                if(xIsInf<IssmDouble>(lambday[i])) _error_("Inf found in solution vector");
    23472351                if(dim==3){
    2348                         lambdaz[i] = values[i*dim+2]; if(xIsNan<IssmDouble>(lambdaz[i])) _error_("NaN found in solution vector");
     2352                        lambdaz[i] = values[i*dim+2];
     2353                        if(xIsNan<IssmDouble>(lambdaz[i])) _error_("NaN found in solution vector");
     2354                        if(xIsInf<IssmDouble>(lambdaz[i])) _error_("Inf found in solution vector");
    23492355                }
    23502356        }
    23512357        for(i=0;i<pnumnodes;i++){
    2352                 lambdap[i] = values[vnumdof+i]; if(xIsNan<IssmDouble>(lambdap[i])) _error_("NaN found in solution vector");
     2358                lambdap[i] = values[vnumdof+i];
     2359                if(xIsNan<IssmDouble>(lambdap[i])) _error_("NaN found in solution vector");
     2360                if(xIsInf<IssmDouble>(lambdap[i])) _error_("Inf found in solution vector");
    23532361        }
    23542362
    … …  
    24092417                /*Check solution*/
    24102418                if(xIsNan<IssmDouble>(lambdax[i])) _error_("NaN found in solution vector");
     2419                if(xIsInf<IssmDouble>(lambdax[i])) _error_("Inf found in solution vector");
    24112420                if(domaintype!=Domain2DverticalEnum && xIsNan<IssmDouble>(lambday[i])) _error_("NaN found in solution vector");
     2421                if(domaintype!=Domain2DverticalEnum && xIsInf<IssmDouble>(lambday[i])) _error_("Inf found in solution vector");
    24122422        }
    24132423

  • issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.cpp

    r20661 r20669  
    649649                /*Check solution*/
    650650                if(xIsNan<IssmDouble>(newdamage[i])) _error_("NaN found in solution vector");
     651                if(xIsInf<IssmDouble>(newdamage[i])) _error_("Inf found in solution vector");
    651652                /*Enforce D < max_damage and D > 0 */
    652653                if(newdamage[i]>max_damage) newdamage[i]=max_damage;

  • issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp

    r20645 r20669  
    14661466                /*Check solution*/
    14671467                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     1468                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    14681469        }
    14691470

  • issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp

    r20645 r20669  
    398398                eplHeads[i]=solution[doflist[i]];
    399399                if(xIsNan<IssmDouble>(eplHeads[i])) _error_("NaN found in solution vector");
     400                if(xIsInf<IssmDouble>(eplHeads[i])) _error_("Inf found in solution vector");
    400401        }
    401402        /*Add input to the element: */

  • issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.cpp

    r20645 r20669  
    435435                values[i] =solution[doflist[i]];
    436436                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     437                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    437438        }
    438439

  • issm/trunk-jpl/src/c/analyses/HydrologyShreveAnalysis.cpp

    r20645 r20669  
    356356                values[i]=solution[doflist[i]];
    357357                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     358                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    358359                if (values[i]<10e-10) values[i]=10e-10; //correcting the water column to positive values
    359360        }

  • issm/trunk-jpl/src/c/analyses/HydrologySommersAnalysis.cpp

    r20645 r20669  
    345345
    346346                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     347                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    347348        }
    348349

  • issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp

    r20655 r20669  
    684684                newthickness[i]=solution[doflist[i]];
    685685                if(xIsNan<IssmDouble>(newthickness[i])) _error_("NaN found in solution vector");
     686                if(xIsInf<IssmDouble>(newthickness[i])) _error_("Inf found in solution vector");
    686687                /*Constrain thickness to be at least 1m*/
    687688                if(newthickness[i]<minthickness) newthickness[i]=minthickness;

  • issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp

    r20645 r20669  
    18511851                vx[i]=values[i*dim+0];
    18521852                if(xIsNan<IssmDouble>(vx[i])) _error_("NaN found in solution vector");
     1853                if(xIsInf<IssmDouble>(vx[i])) _error_("Inf found in solution vector");
    18531854
    18541855                if(dim==2){
    18551856                        vy[i]=values[i*dim+1];
    18561857                        if(xIsNan<IssmDouble>(vy[i])) _error_("NaN found in solution vector");
     1858                        if(xIsInf<IssmDouble>(vy[i])) _error_("Inf found in solution vector");
    18571859                }
    18581860        }
    … …  
    21872189                /*Check solution*/
    21882190                if(xIsNan<IssmDouble>(vx[i])) _error_("NaN found in solution vector");
     2191                if(xIsInf<IssmDouble>(vx[i])) _error_("Inf found in solution vector");
    21892192                if(xIsNan<IssmDouble>(vy[i])) _error_("NaN found in solution vector");
     2193                if(xIsInf<IssmDouble>(vy[i])) _error_("Inf found in solution vector");
    21902194        }
    21912195
    … …  
    28152819                vx[i]=values[i*(dim-1)+0];
    28162820                if(xIsNan<IssmDouble>(vx[i])) _error_("NaN found in solution vector");
     2821                if(xIsInf<IssmDouble>(vx[i])) _error_("Inf found in solution vector");
    28172822                if(dim==3){
    28182823                        vy[i]=values[i*(dim-1)+1];
    28192824                        if(xIsNan<IssmDouble>(vy[i])) _error_("NaN found in solution vector");
     2825                        if(xIsInf<IssmDouble>(vy[i])) _error_("Inf found in solution vector");
    28202826                }
    28212827        }
    … …  
    49534959                vy[i] = values[i*dim+1];
    49544960                if(xIsNan<IssmDouble>(vx[i])) _error_("NaN found in solution vector");
     4961                if(xIsInf<IssmDouble>(vx[i])) _error_("Inf found in solution vector");
    49554962                if(xIsNan<IssmDouble>(vy[i])) _error_("NaN found in solution vector");
     4963                if(xIsInf<IssmDouble>(vy[i])) _error_("Inf found in solution vector");
    49564964
    49574965                if(dim==3){
    49584966                        vz[i] = values[i*dim+2];
    49594967                        if(xIsNan<IssmDouble>(vz[i])) _error_("NaN found in solution vector");
     4968                        if(xIsInf<IssmDouble>(vz[i])) _error_("Inf found in solution vector");
    49604969                }
    49614970        }
    … …  
    49634972                pressure[i] = values[vnumdof+i];
    49644973                if(xIsNan<IssmDouble>(pressure[i])) _error_("NaN found in solution vector");
     4974                if(xIsInf<IssmDouble>(pressure[i])) _error_("Inf found in solution vector");
    49654975        }
    49664976
    … …  
    68546864                /*Check solution*/
    68556865                if(xIsNan<IssmDouble>(vx[i]))       _error_("NaN found in solution vector");
     6866                if(xIsInf<IssmDouble>(vx[i]))       _error_("Inf found in solution vector");
    68566867                if(xIsNan<IssmDouble>(vy[i]))       _error_("NaN found in solution vector");
     6868                if(xIsInf<IssmDouble>(vy[i]))       _error_("Inf found in solution vector");
    68576869                if(xIsNan<IssmDouble>(vzFS[i]))     _error_("NaN found in solution vector");
     6870                if(xIsInf<IssmDouble>(vzFS[i]))     _error_("Inf found in solution vector");
    68586871                if(xIsNan<IssmDouble>(pressure[i])) _error_("NaN found in solution vector");
     6872                if(xIsInf<IssmDouble>(pressure[i])) _error_("Inf found in solution vector");
    68596873        }
    68606874
    … …  
    69606974                /*Check solution*/
    69616975                if(xIsNan<IssmDouble>(vx[i]))       _error_("NaN found in solution vector");
     6976                if(xIsInf<IssmDouble>(vx[i]))       _error_("Inf found in solution vector");
    69626977                if(xIsNan<IssmDouble>(vy[i]))       _error_("NaN found in solution vector");
     6978                if(xIsInf<IssmDouble>(vy[i]))       _error_("Inf found in solution vector");
    69636979                if(xIsNan<IssmDouble>(vzFS[i]))     _error_("NaN found in solution vector");
     6980                if(xIsInf<IssmDouble>(vzFS[i]))     _error_("Inf found in solution vector");
    69646981                if(xIsNan<IssmDouble>(pressure[i])) _error_("NaN found in solution vector");
     6982                if(xIsInf<IssmDouble>(pressure[i])) _error_("Inf found in solution vector");
    69656983        }
    69666984
    … …  
    70557073                /*Check solution*/
    70567074                if(xIsNan<IssmDouble>(vx[i])) _error_("NaN found in solution vector");
     7075                if(xIsInf<IssmDouble>(vx[i])) _error_("Inf found in solution vector");
    70577076                if(xIsNan<IssmDouble>(vy[i])) _error_("NaN found in solution vector");
     7077                if(xIsInf<IssmDouble>(vy[i])) _error_("Inf found in solution vector");
    70587078        }
    70597079

  • issm/trunk-jpl/src/c/analyses/StressbalanceSIAAnalysis.cpp

    r20645 r20669  
    586586                /*Check solution*/
    587587                if(xIsNan<IssmDouble>(vx[i])) _error_("NaN found in solution vector");
     588                if(xIsInf<IssmDouble>(vx[i])) _error_("Inf found in solution vector");
    588589                if(xIsNan<IssmDouble>(vy[i])) _error_("NaN found in solution vector");
     590                if(xIsInf<IssmDouble>(vy[i])) _error_("Inf found in solution vector");
    589591        }
    590592

  • issm/trunk-jpl/src/c/analyses/StressbalanceVerticalAnalysis.cpp

    r20645 r20669  
    539539                /*Check solution*/
    540540                if(xIsNan<IssmDouble>(vz[i])) _error_("NaN found in solution vector");
     541                if(xIsInf<IssmDouble>(vz[i])) _error_("Inf found in solution vector");
    541542        }
    542543

  • issm/trunk-jpl/src/c/analyses/ThermalAnalysis.cpp

    r20645 r20669  
    739739                /*Check solution*/
    740740                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     741                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    741742                //if(values[i]<0)      _printf_("temperature < 0°K found in solution vector\n");
    742743                //if(values[i]>275)    _printf_("temperature > 275°K found in solution vector (Paterson's rheology associated is negative)\n");

  • issm/trunk-jpl/src/c/classes/Elements/Penta.cpp

    r20663 r20669  
    15061506                values[i]=solution[doflist[i]];
    15071507                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     1508                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    15081509        }
    15091510
    … …  
    15361537                values[i+numdof2d]=values[i];
    15371538                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     1539                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    15381540        }
    15391541
    … …  
    15891591                                values[i]=vector[doflist[i]];
    15901592                                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     1593                                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    15911594                        }
    15921595                        /*Add input to the element: */
    … …  
    16011604                                values[i]=vector[nodes[i]->Sid()];
    16021605                                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     1606                                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    16031607                        }
    16041608                        /*Add input to the element: */

  • issm/trunk-jpl/src/c/classes/Elements/Tetra.cpp

    r20645 r20669  
    475475                values[i]=solution[doflist[i]];
    476476                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     477                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    477478        }
    478479

  • issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

    r20659 r20669  
    19271927                values[i]=solution[doflist[i]];
    19281928                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
     1929                if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
    19291930        }
    19301931
    … …  
    19741975                        values[i]=vector[doflist[i]];
    19751976                        if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in vector");
     1977                        if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in vector");
    19761978                }
    19771979                this->inputs->AddInput(new TriaInput(name,values,this->element_type));
    … …  
    19861988                        values[i]=vector[nodes[i]->Sid()];
    19871989                        if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in vector");
     1990                        if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in vector");
    19881991                }
    19891992                this->inputs->AddInput(new TriaInput(name,values,this->element_type));

  • issm/trunk-jpl/src/c/classes/matrix/ElementMatrix.cpp

    r18344 r20669  
    368368                for(int j=0;j<this->ncols;j++){
    369369                        if (xIsNan<IssmDouble>(this->values[i*this->ncols+j])) _error_("NaN found in Element Matrix");
     370                        if (xIsInf<IssmDouble>(this->values[i*this->ncols+j])) _error_("Inf found in Element Matrix");
    370371                        if (fabs(this->values[i*this->ncols+j])>1.e+50) _error_("Element Matrix values exceeds 1.e+50");
    371372                }

  • issm/trunk-jpl/src/c/classes/matrix/ElementVector.cpp

    r18063 r20669  
    201201        for (int i=0;i<this->nrows;i++){
    202202                if (xIsNan<IssmDouble>(this->values[i])) _error_("NaN found in Element Vector");
     203                if (xIsInf<IssmDouble>(this->values[i])) _error_("Inf found in Element Vector");
    203204                if (fabs( this->values[i])>1.e+50) _error_("Element Vector values exceeds 1.e+50");
    204205        }

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

    r18125 r20669  
    77#endif
    88
    9 #include "isnan.h"
     9#include "./isnan.h"
    1010
    1111#if defined(_HAVE_ADOLC_) && !defined(_WRAPPERS_)
    … …  
    1414}
    1515#endif
     16
     17#if defined(_HAVE_ADOLC_) && !defined(_WRAPPERS_)
     18template <> int xIsInf<adouble> (const adouble& X){
     19  return isinf(X.getValue());
     20}
     21#endif

  • issm/trunk-jpl/src/c/shared/Numerics/isnan.h

    r19490 r20669  
    2323}
    2424
     25template <class T> int xIsInf(const T& X) {
     26        return std::isinf(X);
     27}
     28
    2529#if defined(_HAVE_ADOLC_) && !defined(_WRAPPERS_)
    2630#include "./types.h"
    2731template <> int xIsNan<adouble> (const adouble& X);
     32template <> int xIsInf<adouble> (const adouble& X);
    2833#endif
    2934
Note: See TracChangeset for help on using the changeset viewer.