Changeset 12472


Ignore:
Timestamp:
06/20/12 09:19:54 (13 years ago)
Author:
utke
Message:

type renames

Location:
issm/trunk-jpl/src/c/objects
Files:
24 edited

Legend:

Unmodified
Added
Removed

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

    r12455 r12472  
    3434        this->inputs=inputs_in;
    3535        this->element_type=xNew<char>(strlen(element_type_in)+1);
    36         memcpy(this->element_type,element_type_in,(strlen(element_type_in)+1)*sizeof(char));
     36        xMemCpy<char>(this->element_type,element_type_in,(strlen(element_type_in)+1));
    3737
    3838        this->matpar=matpar_in;

  • issm/trunk-jpl/src/c/objects/Loads/Icefront.cpp

    r12455 r12472  
    279279/*}}}*/
    280280/*FUNCTION Icefront::PenaltyCreateKMatrix {{{*/
    281 void  Icefront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, double kmax){
     281void  Icefront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, IssmDouble kmax){
    282282        /*do nothing: */
    283283        return;
    … …  
    285285/*}}}*/
    286286/*FUNCTION Icefront::PenaltyCreatePVector{{{*/
    287 void  Icefront::PenaltyCreatePVector(Vector* pf,double kmax){
     287void  Icefront::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){
    288288        /*do nothing: */
    289289        return;
    … …  
    291291/*}}}*/
    292292/*FUNCTION Icefront::PenaltyCreateJacobianMatrix{{{*/
    293 void  Icefront::PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){
     293void  Icefront::PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){
    294294        this->PenaltyCreateKMatrix(Jff,NULL,kmax);
    295295}
    … …  
    303303
    304304/*Update virtual functions definitions:*/
    305 /*FUNCTION Icefront::InputUpdateFromVector(double* vector, int name, int type) {{{*/
    306 void  Icefront::InputUpdateFromVector(double* vector, int name, int type){
     305/*FUNCTION Icefront::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/
     306void  Icefront::InputUpdateFromVector(IssmDouble* vector, int name, int type){
    307307        /*Nothing updated yet*/
    308308}
    … …  
    318318}
    319319/*}}}*/
    320 /*FUNCTION Icefront::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols, int name, int type) {{{*/
    321 void  Icefront::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols, int name, int type){
    322         /*Nothing updated yet*/
    323 }
    324 /*}}}*/
    325 /*FUNCTION Icefront::InputUpdateFromVectorDakota(double* vector, int name, int type) {{{*/
    326 void  Icefront::InputUpdateFromVectorDakota(double* vector, int name, int type){
     320/*FUNCTION Icefront::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type) {{{*/
     321void  Icefront::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type){
     322        /*Nothing updated yet*/
     323}
     324/*}}}*/
     325/*FUNCTION Icefront::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type) {{{*/
     326void  Icefront::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){
    327327        /*Nothing updated yet*/
    328328}
    … …  
    338338}
    339339/*}}}*/
    340 /*FUNCTION Icefront::InputUpdateFromConstant(double constant, int name) {{{*/
    341 void  Icefront::InputUpdateFromConstant(double constant, int name){
     340/*FUNCTION Icefront::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/
     341void  Icefront::InputUpdateFromConstant(IssmDouble constant, int name){
    342342        /*Nothing updated yet*/
    343343}
    … …  
    354354/*}}}*/
    355355/*FUNCTION Icefront::InputUpdateFromSolution{{{*/
    356 void  Icefront::InputUpdateFromSolution(double* solution){
     356void  Icefront::InputUpdateFromSolution(IssmDouble* solution){
    357357        /*Nothing updated yet*/
    358358}
    … …  
    392392        /*Intermediary*/
    393393        int        ig,index1,index2,fill;
    394         double     Jdet;
    395         double     thickness,bed,pressure,ice_pressure,rho_water,rho_ice,gravity;
    396         double     water_pressure,air_pressure,surface_under_water,base_under_water;
    397         double     xyz_list[numnodes][3];
    398         double     normal[2];
    399         double     L[2];
     394        IssmDouble     Jdet;
     395        IssmDouble     thickness,bed,pressure,ice_pressure,rho_water,rho_ice,gravity;
     396        IssmDouble     water_pressure,air_pressure,surface_under_water,base_under_water;
     397        IssmDouble     xyz_list[numnodes][3];
     398        IssmDouble     normal[2];
     399        IssmDouble     L[2];
    400400        GaussTria *gauss;
    401401
    … …  
    513513        int         i,j,ig,index1,index2,index3,index4;
    514514        int         fill;
    515         double      surface,pressure,ice_pressure,rho_water,rho_ice,gravity;
    516         double      water_pressure,air_pressure;
    517         double      Jdet,z_g;
    518         double      xyz_list[NUMVERTICESQUA][3];
    519         double      normal[3];
    520         double      l1l4[4];
     515        IssmDouble      surface,pressure,ice_pressure,rho_water,rho_ice,gravity;
     516        IssmDouble      water_pressure,air_pressure;
     517        IssmDouble      Jdet,z_g;
     518        IssmDouble      xyz_list[NUMVERTICESQUA][3];
     519        IssmDouble      normal[3];
     520        IssmDouble      l1l4[4];
    521521        GaussPenta *gauss = NULL;
    522522
    … …  
    543543
    544544        /* Start  looping on the number of gaussian points: */
    545         double zmax=xyz_list[0][2]; for(i=1;i<NUMVERTICESQUA;i++) if(xyz_list[i][2]>zmax) zmax=xyz_list[i][2];
    546         double zmin=xyz_list[0][2]; for(i=1;i<NUMVERTICESQUA;i++) if(xyz_list[i][2]<zmin) zmin=xyz_list[i][2];
     545        IssmDouble zmax=xyz_list[0][2]; for(i=1;i<NUMVERTICESQUA;i++) if(xyz_list[i][2]>zmax) zmax=xyz_list[i][2];
     546        IssmDouble zmin=xyz_list[0][2]; for(i=1;i<NUMVERTICESQUA;i++) if(xyz_list[i][2]<zmin) zmin=xyz_list[i][2];
    547547        if(zmax>0 && zmin<0) gauss=new GaussPenta(index1,index2,index3,index4,3,10); //refined in vertical because of the sea level discontinuity
    548548        else                 gauss=new GaussPenta(index1,index2,index3,index4,3,3);
    … …  
    590590        int         i,j,ig,index1,index2,index3,index4;
    591591        int         fill;
    592         double      pressure,rho_water,gravity;
    593         double      water_pressure,air_pressure;
    594         double      Jdet,z_g;
    595         double      xyz_list[NUMVERTICESQUA][3];
    596         double      normal[3];
    597         double      l1l4[4];
     592        IssmDouble      pressure,rho_water,gravity;
     593        IssmDouble      water_pressure,air_pressure;
     594        IssmDouble      Jdet,z_g;
     595        IssmDouble      xyz_list[NUMVERTICESQUA][3];
     596        IssmDouble      normal[3];
     597        IssmDouble      l1l4[4];
    598598        GaussPenta *gauss = NULL;
    599599
    … …  
    618618
    619619        /* Start  looping on the number of gaussian points: */
    620         double zmax=xyz_list[0][2]; for(i=1;i<NUMVERTICESQUA;i++) if(xyz_list[i][2]>zmax) zmax=xyz_list[i][2];
    621         double zmin=xyz_list[0][2]; for(i=1;i<NUMVERTICESQUA;i++) if(xyz_list[i][2]<zmin) zmin=xyz_list[i][2];
     620        IssmDouble zmax=xyz_list[0][2]; for(i=1;i<NUMVERTICESQUA;i++) if(xyz_list[i][2]>zmax) zmax=xyz_list[i][2];
     621        IssmDouble zmin=xyz_list[0][2]; for(i=1;i<NUMVERTICESQUA;i++) if(xyz_list[i][2]<zmin) zmin=xyz_list[i][2];
    622622        if(zmax>0 && zmin<0) gauss=new GaussPenta(index1,index2,index3,index4,3,30); //refined in vertical because of the sea level discontinuity
    623623        else                 gauss=new GaussPenta(index1,index2,index3,index4,3,3);
    … …  
    705705/*}}}*/
    706706/*FUNCTION Icefront::GetSegmentNormal {{{*/
    707 void Icefront:: GetSegmentNormal(double* normal,double xyz_list[4][3]){
     707void Icefront:: GetSegmentNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]){
    708708
    709709        /*Build unit outward pointing vector*/
    710710        const int numnodes=NUMVERTICESSEG;
    711         double vector[2];
    712         double norm;
     711        IssmDouble vector[2];
     712        IssmDouble norm;
    713713
    714714        vector[0]=xyz_list[1][0] - xyz_list[0][0];
    … …  
    722722/*}}}*/
    723723/*FUNCTION Icefront::GetQuadNormal {{{*/
    724 void Icefront:: GetQuadNormal(double* normal,double xyz_list[4][3]){
     724void Icefront:: GetQuadNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]){
    725725
    726726        /*Build unit outward pointing vector*/
    727         double AB[3];
    728         double AC[3];
    729         double norm;
     727        IssmDouble AB[3];
     728        IssmDouble AC[3];
     729        IssmDouble norm;
    730730
    731731        AB[0]=xyz_list[1][0] - xyz_list[0][0];

  • issm/trunk-jpl/src/c/objects/Loads/Icefront.h

    r12365 r12472  
    5353                /*}}}*/
    5454                /*Update virtual functions definitions: {{{*/
    55                 void  InputUpdateFromVector(double* vector, int name, int type);
     55                void  InputUpdateFromVector(IssmDouble* vector, int name, int type);
    5656                void  InputUpdateFromVector(int* vector, int name, int type);
    5757                void  InputUpdateFromVector(bool* vector, int name, int type);
    58                 void  InputUpdateFromMatrixDakota(double* matrix,int ncols,int nrows, int name, int type);
    59                 void  InputUpdateFromVectorDakota(double* vector, int name, int type);
     58                void  InputUpdateFromMatrixDakota(IssmDouble* matrix,int ncols,int nrows, int name, int type);
     59                void  InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type);
    6060                void  InputUpdateFromVectorDakota(int* vector, int name, int type);
    6161                void  InputUpdateFromVectorDakota(bool* vector, int name, int type);
    62                 void  InputUpdateFromConstant(double constant, int name);
     62                void  InputUpdateFromConstant(IssmDouble constant, int name);
    6363                void  InputUpdateFromConstant(int constant, int name);
    6464                void  InputUpdateFromConstant(bool constant, int name);
    65                 void  InputUpdateFromSolution(double* solution);
     65                void  InputUpdateFromSolution(IssmDouble* solution);
    6666                void  InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");};
    6767                /*}}}*/
    … …  
    7272                void  CreatePVector(Vector* pf);
    7373                void  CreateJacobianMatrix(Matrix* Jff);
    74                 void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax);
    75                 void  PenaltyCreatePVector(Vector*  pf, double kmax);
    76                 void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax);
     74                void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, IssmDouble kmax);
     75                void  PenaltyCreatePVector(Vector*  pf, IssmDouble kmax);
     76                void  PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax);
    7777                bool  InAnalysis(int analysis_type);
    7878                /*}}}*/
    7979                /*Load management: {{{*/
    8080                void GetDofList(int** pdoflist,int approximation_enum,int setenum);
    81                 void GetSegmentNormal(double* normal,double xyz_list[2][3]);
    82                 void GetQuadNormal(double* normal,double xyz_list[4][3]);
     81                void GetSegmentNormal(IssmDouble* normal,IssmDouble xyz_list[2][3]);
     82                void GetQuadNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]);
    8383                #ifdef _HAVE_CONTROL_
    8484                ElementVector* CreatePVectorAdjointHoriz(void);

  • issm/trunk-jpl/src/c/objects/Loads/Load.h

    r12365 r12472  
    3232                virtual void  CreatePVector(Vector* pf)=0;
    3333                virtual void  CreateJacobianMatrix(Matrix* Jff)=0;
    34                 virtual void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax)=0;
    35                 virtual void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, double kmax)=0;
    36                 virtual void  PenaltyCreatePVector(Vector* pf, double kmax)=0;
     34                virtual void  PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax)=0;
     35                virtual void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, IssmDouble kmax)=0;
     36                virtual void  PenaltyCreatePVector(Vector* pf, IssmDouble kmax)=0;
    3737                virtual bool  InAnalysis(int analysis_type)=0;
    3838                /*}}}*/

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

    r12365 r12472  
    6060
    6161        /*First, see wether this is an internal or boundary edge (if e2=NaN)*/
    62         if (isnan((double)iomodel->Data(MeshEdgesEnum)[4*i+3])){ //edges are [node1 node2 elem1 elem2]
     62        if (isnan((IssmDouble)iomodel->Data(MeshEdgesEnum)[4*i+3])){ //edges are [node1 node2 elem1 elem2]
    6363                /* Boundary edge, only one element */
    6464                e1=(int)iomodel->Data(MeshEdgesEnum)[4*i+2];
    … …  
    312312/*}}}*/
    313313/*FUNCTION Numericalflux::PenaltyCreateKMatrix {{{*/
    314 void  Numericalflux::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){
     314void  Numericalflux::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,IssmDouble kmax){
    315315
    316316        /*No stiffness loads applied, do nothing: */
    … …  
    320320/*}}}*/
    321321/*FUNCTION Numericalflux::PenaltyCreatePVector{{{*/
    322 void  Numericalflux::PenaltyCreatePVector(Vector* pf,double kmax){
     322void  Numericalflux::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){
    323323
    324324        /*No penalty loads applied, do nothing: */
    … …  
    359359        /* Intermediaries*/
    360360        int        i,j,ig,index1,index2;
    361         double     DL1,DL2,Jdet,dt,vx,vy,UdotN;
    362         double     xyz_list[NUMVERTICES_INTERNAL][3];
    363         double     normal[2];
    364         double     B[numdof];
    365         double     Bprime[numdof];
    366         double     Ke_g1[numdof][numdof];
    367         double     Ke_g2[numdof][numdof];
     361        IssmDouble     DL1,DL2,Jdet,dt,vx,vy,UdotN;
     362        IssmDouble     xyz_list[NUMVERTICES_INTERNAL][3];
     363        IssmDouble     normal[2];
     364        IssmDouble     B[numdof];
     365        IssmDouble     Bprime[numdof];
     366        IssmDouble     Ke_g1[numdof][numdof];
     367        IssmDouble     Ke_g2[numdof][numdof];
    368368        GaussTria *gauss;
    369369
    … …  
    424424        /* Intermediaries*/
    425425        int        i,j,ig,index1,index2;
    426         double     DL,Jdet,dt,vx,vy,mean_vx,mean_vy,UdotN;
    427         double     xyz_list[NUMVERTICES_BOUNDARY][3];
    428         double     normal[2];
    429         double     L[numdof];
    430         double     Ke_g[numdof][numdof];
     426        IssmDouble     DL,Jdet,dt,vx,vy,mean_vx,mean_vy,UdotN;
     427        IssmDouble     xyz_list[NUMVERTICES_BOUNDARY][3];
     428        IssmDouble     normal[2];
     429        IssmDouble     L[numdof];
     430        IssmDouble     Ke_g[numdof][numdof];
    431431        GaussTria *gauss;
    432432
    … …  
    512512        /* Intermediaries*/
    513513        int        i,j,ig,index1,index2;
    514         double     DL1,DL2,Jdet,vx,vy,UdotN;
    515         double     xyz_list[NUMVERTICES_INTERNAL][3];
    516         double     normal[2];
    517         double     B[numdof];
    518         double     Bprime[numdof];
    519         double     Ke_g1[numdof][numdof];
    520         double     Ke_g2[numdof][numdof];
     514        IssmDouble     DL1,DL2,Jdet,vx,vy,UdotN;
     515        IssmDouble     xyz_list[NUMVERTICES_INTERNAL][3];
     516        IssmDouble     normal[2];
     517        IssmDouble     B[numdof];
     518        IssmDouble     Bprime[numdof];
     519        IssmDouble     Ke_g1[numdof][numdof];
     520        IssmDouble     Ke_g2[numdof][numdof];
    521521        GaussTria *gauss;
    522522
    … …  
    576576        /* Intermediaries*/
    577577        int        i,j,ig,index1,index2;
    578         double     DL,Jdet,vx,vy,mean_vx,mean_vy,UdotN;
    579         double     xyz_list[NUMVERTICES_BOUNDARY][3];
    580         double     normal[2];
    581         double     L[numdof];
    582         double     Ke_g[numdof][numdof];
     578        IssmDouble     DL,Jdet,vx,vy,mean_vx,mean_vy,UdotN;
     579        IssmDouble     xyz_list[NUMVERTICES_BOUNDARY][3];
     580        IssmDouble     normal[2];
     581        IssmDouble     L[numdof];
     582        IssmDouble     Ke_g[numdof][numdof];
    583583        GaussTria *gauss;
    584584
    … …  
    703703        /* Intermediaries*/
    704704        int        i,j,ig,index1,index2;
    705         double     DL,Jdet,dt,vx,vy,mean_vx,mean_vy,UdotN,thickness;
    706         double     xyz_list[NUMVERTICES_BOUNDARY][3];
    707         double     normal[2];
    708         double     L[numdof];
     705        IssmDouble     DL,Jdet,dt,vx,vy,mean_vx,mean_vy,UdotN,thickness;
     706        IssmDouble     xyz_list[NUMVERTICES_BOUNDARY][3];
     707        IssmDouble     normal[2];
     708        IssmDouble     L[numdof];
    709709        GaussTria *gauss;
    710710
    … …  
    797797        /* Intermediaries*/
    798798        int        i,j,ig,index1,index2;
    799         double     DL,Jdet,vx,vy,mean_vx,mean_vy,UdotN,thickness;
    800         double     xyz_list[NUMVERTICES_BOUNDARY][3];
    801         double     normal[2];
    802         double     L[numdof];
     799        IssmDouble     DL,Jdet,vx,vy,mean_vx,mean_vy,UdotN,thickness;
     800        IssmDouble     xyz_list[NUMVERTICES_BOUNDARY][3];
     801        IssmDouble     normal[2];
     802        IssmDouble     L[numdof];
    803803        GaussTria *gauss;
    804804
    … …  
    864864/*}}}*/
    865865/*FUNCTION Numericalflux::GetNormal {{{*/
    866 void Numericalflux:: GetNormal(double* normal,double xyz_list[4][3]){
     866void Numericalflux:: GetNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]){
    867867
    868868        /*Build unit outward pointing vector*/
    869         double vector[2];
    870         double norm;
     869        IssmDouble vector[2];
     870        IssmDouble norm;
    871871
    872872        vector[0]=xyz_list[1][0] - xyz_list[0][0];

  • issm/trunk-jpl/src/c/objects/Loads/Numericalflux.h

    r12365 r12472  
    4949                /*}}}*/
    5050                /*Update virtual functions resolution: {{{*/
    51                 void    InputUpdateFromVector(double* vector, int name, int type){/*Do nothing*/}
     51                void    InputUpdateFromVector(IssmDouble* vector, int name, int type){/*Do nothing*/}
    5252                void    InputUpdateFromVector(int* vector, int name, int type){_error_("Not implemented yet!");}
    5353                void    InputUpdateFromVector(bool* vector, int name, int type){_error_("Not implemented yet!");}
    54                 void    InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols, int name, int type){/*Do nothing*/}
    55                 void    InputUpdateFromVectorDakota(double* vector, int name, int type){/*Do nothing*/}
     54                void    InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type){/*Do nothing*/}
     55                void    InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){/*Do nothing*/}
    5656                void    InputUpdateFromVectorDakota(int* vector, int name, int type){_error_("Not implemented yet!");}
    5757                void    InputUpdateFromVectorDakota(bool* vector, int name, int type){_error_("Not implemented yet!");}
    58                 void    InputUpdateFromConstant(double constant, int name){/*Do nothing*/};
     58                void    InputUpdateFromConstant(IssmDouble constant, int name){/*Do nothing*/};
    5959                void    InputUpdateFromConstant(int constant, int name){/*Do nothing*/};
    6060                void    InputUpdateFromConstant(bool constant, int name){_error_("Not implemented yet!");}
    61                 void    InputUpdateFromSolution(double* solution){_error_("Not implemented yet!");}
     61                void    InputUpdateFromSolution(IssmDouble* solution){_error_("Not implemented yet!");}
    6262                void  InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");};
    6363                /*}}}*/
    … …  
    6868                void  CreatePVector(Vector* pf);
    6969                void  CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");};
    70                 void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");};
    71                 void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax);
    72                 void  PenaltyCreatePVector(Vector* pf, double kmax);
     70                void  PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){_error_("Not implemented yet");};
     71                void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, IssmDouble kmax);
     72                void  PenaltyCreatePVector(Vector* pf, IssmDouble kmax);
    7373                bool  InAnalysis(int analysis_type);
    7474                /*}}}*/
    7575                /*Numericalflux management:{{{*/
    76                 void  GetNormal(double* normal,double xyz_list[4][3]);
     76                void  GetNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]);
    7777                ElementMatrix* CreateKMatrixPrognostic(void);
    7878                ElementMatrix* CreateKMatrixPrognosticInternal(void);

  • issm/trunk-jpl/src/c/objects/Loads/Pengrid.cpp

    r12365 r12472  
    217217/*}}}*/
    218218/*FUNCTION Pengrid::PenaltyCreateMatrix {{{*/
    219 void  Pengrid::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){
     219void  Pengrid::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,IssmDouble kmax){
    220220
    221221        /*Retrieve parameters: */
    … …  
    250250/*}}}*/
    251251/*FUNCTION Pengrid::PenaltyCreatePVector {{{*/
    252 void  Pengrid::PenaltyCreatePVector(Vector* pf,double kmax){
     252void  Pengrid::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){
    253253
    254254        /*Retrieve parameters: */
    … …  
    289289
    290290/*Update virtual functions definitions:*/
    291 /*FUNCTION Pengrid::InputUpdateFromVector(double* vector, int name, int type) {{{*/
    292 void  Pengrid::InputUpdateFromVector(double* vector, int name, int type){
     291/*FUNCTION Pengrid::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/
     292void  Pengrid::InputUpdateFromVector(IssmDouble* vector, int name, int type){
    293293        /*Nothing updated yet*/
    294294}
    … …  
    304304}
    305305/*}}}*/
    306 /*FUNCTION Pengrid::InputUpdateFromMatrixDakota(double* vector, int nrows, int ncols, int name, int type) {{{*/
    307 void  Pengrid::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols, int name, int type){
    308         /*Nothing updated yet*/
    309 }
    310 /*}}}*/
    311 /*FUNCTION Pengrid::InputUpdateFromVectorDakota(double* vector, int name, int type) {{{*/
    312 void  Pengrid::InputUpdateFromVectorDakota(double* vector, int name, int type){
     306/*FUNCTION Pengrid::InputUpdateFromMatrixDakota(IssmDouble* vector, int nrows, int ncols, int name, int type) {{{*/
     307void  Pengrid::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type){
     308        /*Nothing updated yet*/
     309}
     310/*}}}*/
     311/*FUNCTION Pengrid::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type) {{{*/
     312void  Pengrid::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){
    313313        /*Nothing updated yet*/
    314314}
    … …  
    324324}
    325325/*}}}*/
    326 /*FUNCTION Pengrid::InputUpdateFromConstant(double constant, int name) {{{*/
    327 void  Pengrid::InputUpdateFromConstant(double constant, int name){
     326/*FUNCTION Pengrid::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/
     327void  Pengrid::InputUpdateFromConstant(IssmDouble constant, int name){
    328328        switch(name){
    329329
    … …  
    353353/*}}}*/
    354354/*FUNCTION Pengrid::InputUpdateFromSolution{{{*/
    355 void  Pengrid::InputUpdateFromSolution(double* solution){
     355void  Pengrid::InputUpdateFromSolution(IssmDouble* solution){
    356356        /*Nothing updated yet*/
    357357}
    … …  
    391391        int    found=0;
    392392        const int numnodes=1;
    393         double pressure;
    394         double temperature;
    395         double t_pmp;
     393        IssmDouble pressure;
     394        IssmDouble temperature;
     395        IssmDouble t_pmp;
    396396        int    new_active;
    397397        int    unstable=0;
    … …  
    455455#ifdef _HAVE_DIAGNOSTIC_
    456456/*FUNCTION Pengrid::PenaltyCreateKMatrixDiagnosticStokes {{{*/
    457 ElementMatrix* Pengrid::PenaltyCreateKMatrixDiagnosticStokes(double kmax){
     457ElementMatrix* Pengrid::PenaltyCreateKMatrixDiagnosticStokes(IssmDouble kmax){
    458458       
    459459        const int numdof = NUMVERTICES *NDOF4;
    460         double    slope[2];
    461         double    penalty_offset;
     460        IssmDouble    slope[2];
     461        IssmDouble    penalty_offset;
    462462        int       approximation;
    463463
    … …  
    475475
    476476        /*Create elementary matrix: add penalty to constrain wb (wb=ub*db/dx+vb*db/dy)*/
    477         Ke->values[2*NDOF4+0]=-slope[0]*kmax*pow((double)10.0,penalty_offset);
    478         Ke->values[2*NDOF4+1]=-slope[1]*kmax*pow((double)10.0,penalty_offset);
    479         Ke->values[2*NDOF4+2]= kmax*pow((double)10,penalty_offset);
     477        Ke->values[2*NDOF4+0]=-slope[0]*kmax*pow((IssmDouble)10.0,penalty_offset);
     478        Ke->values[2*NDOF4+1]=-slope[1]*kmax*pow((IssmDouble)10.0,penalty_offset);
     479        Ke->values[2*NDOF4+2]= kmax*pow((IssmDouble)10,penalty_offset);
    480480
    481481        /*Transform Coordinate System*/
    … …  
    489489#ifdef _HAVE_THERMAL_
    490490/*FUNCTION Pengrid::PenaltyCreateKMatrixMelting {{{*/
    491 ElementMatrix* Pengrid::PenaltyCreateKMatrixMelting(double kmax){
     491ElementMatrix* Pengrid::PenaltyCreateKMatrixMelting(IssmDouble kmax){
    492492
    493493        const int numdof=NUMVERTICES*NDOF1;
    494         double pressure,temperature,t_pmp;
    495         double penalty_factor;
     494        IssmDouble pressure,temperature,t_pmp;
     495        IssmDouble penalty_factor;
    496496
    497497        Penta* penta=(Penta*)element;
    … …  
    511511        /*Add penalty load*/
    512512        if (temperature<t_pmp){ //If T<Tpmp, there must be no melting. Therefore, melting should be  constrained to 0 when T<Tpmp, instead of using spcs, use penalties
    513                 Ke->values[0]=kmax*pow((double)10,penalty_factor);
     513                Ke->values[0]=kmax*pow((IssmDouble)10,penalty_factor);
    514514        }
    515515
    … …  
    519519/*}}}*/
    520520/*FUNCTION Pengrid::PenaltyCreateKMatrixThermal {{{*/
    521 ElementMatrix* Pengrid::PenaltyCreateKMatrixThermal(double kmax){
     521ElementMatrix* Pengrid::PenaltyCreateKMatrixThermal(IssmDouble kmax){
    522522
    523523        const int numdof=NUMVERTICES*NDOF1;
    524         double    penalty_factor;
     524        IssmDouble    penalty_factor;
    525525
    526526        /*Initialize Element matrix and return if necessary*/
    … …  
    531531        parameters->FindParam(&penalty_factor,ThermalPenaltyFactorEnum);
    532532
    533         Ke->values[0]=kmax*pow((double)10,penalty_factor);
     533        Ke->values[0]=kmax*pow((IssmDouble)10,penalty_factor);
    534534
    535535        /*Clean up and return*/
    … …  
    538538/*}}}*/
    539539/*FUNCTION Pengrid::PenaltyCreatePVectorMelting {{{*/
    540 ElementVector* Pengrid::PenaltyCreatePVectorMelting(double kmax){
     540ElementVector* Pengrid::PenaltyCreatePVectorMelting(IssmDouble kmax){
    541541       
    542542        const int numdof=NUMVERTICES*NDOF1;
    543         double pressure;
    544         double temperature;
    545         double melting_offset;
    546         double t_pmp;
    547         double dt,penalty_factor;
     543        IssmDouble pressure;
     544        IssmDouble temperature;
     545        IssmDouble melting_offset;
     546        IssmDouble t_pmp;
     547        IssmDouble dt,penalty_factor;
    548548
    549549        /*recover pointers: */
    … …  
    572572        }
    573573        else{
    574                 if (dt) pe->values[0]=melting_offset*pow((double)10,penalty_factor)*(temperature-t_pmp)/dt;
    575                 else    pe->values[0]=melting_offset*pow((double)10,penalty_factor)*(temperature-t_pmp);
     574                if (dt) pe->values[0]=melting_offset*pow((IssmDouble)10,penalty_factor)*(temperature-t_pmp)/dt;
     575                else    pe->values[0]=melting_offset*pow((IssmDouble)10,penalty_factor)*(temperature-t_pmp);
    576576        }
    577577
    … …  
    581581/*}}}*/
    582582/*FUNCTION Pengrid::PenaltyCreatePVectorThermal {{{*/
    583 ElementVector* Pengrid::PenaltyCreatePVectorThermal(double kmax){
     583ElementVector* Pengrid::PenaltyCreatePVectorThermal(IssmDouble kmax){
    584584
    585585        const int numdof=NUMVERTICES*NDOF1;
    586         double pressure;
    587         double t_pmp;
    588         double penalty_factor;
     586        IssmDouble pressure;
     587        IssmDouble t_pmp;
     588        IssmDouble penalty_factor;
    589589
    590590        Penta* penta=(Penta*)element;
    … …  
    601601        t_pmp=matpar->GetMeltingPoint()-matpar->GetBeta()*pressure;
    602602
    603         pe->values[0]=kmax*pow((double)10,penalty_factor)*t_pmp;
     603        pe->values[0]=kmax*pow((IssmDouble)10,penalty_factor)*t_pmp;
    604604
    605605        /*Clean up and return*/
    … …  
    615615/*}}}*/
    616616/*FUNCTION Pengrid::UpdateInputs {{{*/
    617 void  Pengrid::UpdateInputs(double* solution){
     617void  Pengrid::UpdateInputs(IssmDouble* solution){
    618618        _error_("not supported yet!");
    619619}

  • issm/trunk-jpl/src/c/objects/Loads/Pengrid.h

    r12365 r12472  
    5454                /*}}}*/
    5555                /*Update virtual functions resolution: {{{*/
    56                 void  InputUpdateFromVector(double* vector, int name, int type);
     56                void  InputUpdateFromVector(IssmDouble* vector, int name, int type);
    5757                void  InputUpdateFromVector(int* vector, int name, int type);
    5858                void  InputUpdateFromVector(bool* vector, int name, int type);
    59                 void  InputUpdateFromMatrixDakota(double* matrix ,int nrows, int ncols, int name, int type);
    60                 void  InputUpdateFromVectorDakota(double* vector, int name, int type);
     59                void  InputUpdateFromMatrixDakota(IssmDouble* matrix ,int nrows, int ncols, int name, int type);
     60                void  InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type);
    6161                void  InputUpdateFromVectorDakota(int* vector, int name, int type);
    6262                void  InputUpdateFromVectorDakota(bool* vector, int name, int type);
    63                 void  InputUpdateFromConstant(double constant, int name);
     63                void  InputUpdateFromConstant(IssmDouble constant, int name);
    6464                void  InputUpdateFromConstant(int constant, int name);
    6565                void  InputUpdateFromConstant(bool constant, int name);
    66                 void  InputUpdateFromSolution(double* solution);
     66                void  InputUpdateFromSolution(IssmDouble* solution);
    6767                void  InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");};
    6868                /*}}}*/
    … …  
    7373                void  CreatePVector(Vector* pf);
    7474                void  CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");};
    75                 void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");};
    76                 void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax);
    77                 void  PenaltyCreatePVector(Vector* pf, double kmax);
     75                void  PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){_error_("Not implemented yet");};
     76                void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, IssmDouble kmax);
     77                void  PenaltyCreatePVector(Vector* pf, IssmDouble kmax);
    7878                bool  InAnalysis(int analysis_type);
    7979                /*}}}*/
    8080                /*Pengrid management {{{*/
    8181                #ifdef _HAVE_DIAGNOSTIC_
    82                 ElementMatrix* PenaltyCreateKMatrixDiagnosticStokes(double kmax);
     82                ElementMatrix* PenaltyCreateKMatrixDiagnosticStokes(IssmDouble kmax);
    8383                #endif
    8484                #ifdef _HAVE_THERMAL_
    85                 ElementMatrix* PenaltyCreateKMatrixThermal(double kmax);
    86                 ElementMatrix* PenaltyCreateKMatrixMelting(double kmax);
    87                 ElementVector* PenaltyCreatePVectorThermal(double kmax);
    88                 ElementVector* PenaltyCreatePVectorMelting(double kmax);
     85                ElementMatrix* PenaltyCreateKMatrixThermal(IssmDouble kmax);
     86                ElementMatrix* PenaltyCreateKMatrixMelting(IssmDouble kmax);
     87                ElementVector* PenaltyCreatePVectorThermal(IssmDouble kmax);
     88                ElementVector* PenaltyCreatePVectorMelting(IssmDouble kmax);
    8989                #endif
    9090                void  ConstraintActivate(int* punstable);
    9191                void  ConstraintActivateThermal(int* punstable);
    92                 void  UpdateInputs(double* solution);
     92                void  UpdateInputs(IssmDouble* solution);
    9393                void  ResetConstraint(void);
    9494                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Loads/Penpair.cpp

    r12365 r12472  
    157157/*}}}*/
    158158/*FUNCTION Penpair::PenaltyCreateKMatrix {{{*/
    159 void  Penpair::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){
     159void  Penpair::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,IssmDouble kmax){
    160160
    161161        /*Retrieve parameters: */
    … …  
    183183/*}}}*/
    184184/*FUNCTION Penpair::PenaltyCreatePVector {{{*/
    185 void  Penpair::PenaltyCreatePVector(Vector* pf,double kmax){
     185void  Penpair::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){
    186186        /*No loads applied, do nothing: */
    187187        return;
    … …  
    189189/*}}}*/
    190190/*FUNCTION Penpair::PenaltyCreateJacobianMatrix{{{*/
    191 void  Penpair::PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){
     191void  Penpair::PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){
    192192        this->PenaltyCreateKMatrix(Jff,NULL,kmax);
    193193}
    … …  
    201201
    202202/*Update virtual functions definitions:*/
    203 /*FUNCTION Penpair::InputUpdateFromConstant(double constant, int name) {{{*/
    204 void  Penpair::InputUpdateFromConstant(double constant, int name){
     203/*FUNCTION Penpair::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/
     204void  Penpair::InputUpdateFromConstant(IssmDouble constant, int name){
    205205        /*Nothing updated yet*/
    206206}
    … …  
    216216}
    217217/*}}}*/
    218 /*FUNCTION Penpair::InputUpdateFromVector(double* vector, int name, int type) {{{*/
    219 void  Penpair::InputUpdateFromVector(double* vector, int name, int type){
     218/*FUNCTION Penpair::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/
     219void  Penpair::InputUpdateFromVector(IssmDouble* vector, int name, int type){
    220220        /*Nothing updated yet*/
    221221}
    … …  
    234234/*Penpair management:*/
    235235/*FUNCTION Penpair::PenaltyCreateKMatrixDiagnosticHoriz{{{*/
    236 ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticHoriz(double kmax){
     236ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticHoriz(IssmDouble kmax){
    237237
    238238        int    approximation0=nodes[0]->GetApproximation();
    … …  
    269269/*}}}*/
    270270/*FUNCTION Penpair::PenaltyCreateKMatrixDiagnosticMacAyealPattyn {{{*/
    271 ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticMacAyealPattyn(double kmax){
     271ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticMacAyealPattyn(IssmDouble kmax){
    272272       
    273273        const int numdof=NUMVERTICES*NDOF2;
    274         double penalty_offset;
     274        IssmDouble penalty_offset;
    275275
    276276        /*Initialize Element vector and return if necessary*/
    … …  
    281281
    282282        //Create elementary matrix: add penalty to
    283         Ke->values[0*numdof+0]=+kmax*pow((double)10.0,penalty_offset);
    284         Ke->values[0*numdof+2]=-kmax*pow((double)10.0,penalty_offset);
    285         Ke->values[2*numdof+0]=-kmax*pow((double)10.0,penalty_offset);
    286         Ke->values[2*numdof+2]=+kmax*pow((double)10.0,penalty_offset);
    287 
    288         Ke->values[1*numdof+1]=+kmax*pow((double)10.0,penalty_offset);
    289         Ke->values[1*numdof+3]=-kmax*pow((double)10.0,penalty_offset);
    290         Ke->values[3*numdof+1]=-kmax*pow((double)10.0,penalty_offset);
    291         Ke->values[3*numdof+3]=+kmax*pow((double)10.0,penalty_offset);
     283        Ke->values[0*numdof+0]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     284        Ke->values[0*numdof+2]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     285        Ke->values[2*numdof+0]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     286        Ke->values[2*numdof+2]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     287
     288        Ke->values[1*numdof+1]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     289        Ke->values[1*numdof+3]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     290        Ke->values[3*numdof+1]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     291        Ke->values[3*numdof+3]=+kmax*pow((IssmDouble)10.0,penalty_offset);
    292292
    293293        /*Clean up and return*/
    … …  
    296296/*}}}*/
    297297/*FUNCTION Penpair::PenaltyCreateKMatrixDiagnosticStokes {{{*/
    298 ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticStokes(double kmax){
     298ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticStokes(IssmDouble kmax){
    299299       
    300300        const int numdof=NUMVERTICES*NDOF4;
    301         double penalty_offset;
     301        IssmDouble penalty_offset;
    302302
    303303        /*Initialize Element vector and return if necessary*/
    … …  
    308308
    309309        //Create elementary matrix: add penalty to
    310         Ke->values[0*numdof+0]=+kmax*pow((double)10.0,penalty_offset);
    311         Ke->values[0*numdof+4]=-kmax*pow((double)10.0,penalty_offset);
    312         Ke->values[4*numdof+0]=-kmax*pow((double)10.0,penalty_offset);
    313         Ke->values[4*numdof+4]=+kmax*pow((double)10.0,penalty_offset);
    314 
    315         Ke->values[1*numdof+1]=+kmax*pow((double)10.0,penalty_offset);
    316         Ke->values[1*numdof+5]=-kmax*pow((double)10.0,penalty_offset);
    317         Ke->values[5*numdof+1]=-kmax*pow((double)10.0,penalty_offset);
    318         Ke->values[5*numdof+5]=+kmax*pow((double)10.0,penalty_offset);
    319        
    320         Ke->values[2*numdof+2]=+kmax*pow((double)10.0,penalty_offset);
    321         Ke->values[2*numdof+6]=-kmax*pow((double)10.0,penalty_offset);
    322         Ke->values[6*numdof+2]=-kmax*pow((double)10.0,penalty_offset);
    323         Ke->values[6*numdof+6]=+kmax*pow((double)10.0,penalty_offset);
    324 
    325         Ke->values[3*numdof+3]=+kmax*pow((double)10.0,penalty_offset);
    326         Ke->values[3*numdof+7]=-kmax*pow((double)10.0,penalty_offset);
    327         Ke->values[7*numdof+3]=-kmax*pow((double)10.0,penalty_offset);
    328         Ke->values[7*numdof+7]=+kmax*pow((double)10.0,penalty_offset);
     310        Ke->values[0*numdof+0]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     311        Ke->values[0*numdof+4]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     312        Ke->values[4*numdof+0]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     313        Ke->values[4*numdof+4]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     314
     315        Ke->values[1*numdof+1]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     316        Ke->values[1*numdof+5]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     317        Ke->values[5*numdof+1]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     318        Ke->values[5*numdof+5]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     319       
     320        Ke->values[2*numdof+2]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     321        Ke->values[2*numdof+6]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     322        Ke->values[6*numdof+2]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     323        Ke->values[6*numdof+6]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     324
     325        Ke->values[3*numdof+3]=+kmax*pow((IssmDouble)10.0,penalty_offset);
     326        Ke->values[3*numdof+7]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     327        Ke->values[7*numdof+3]=-kmax*pow((IssmDouble)10.0,penalty_offset);
     328        Ke->values[7*numdof+7]=+kmax*pow((IssmDouble)10.0,penalty_offset);
    329329
    330330        /*Clean up and return*/
    … …  
    333333/*}}}*/
    334334/*FUNCTION Penpair::PenaltyCreateKMatrixPrognostic {{{*/
    335 ElementMatrix* Penpair::PenaltyCreateKMatrixPrognostic(double kmax){
     335ElementMatrix* Penpair::PenaltyCreateKMatrixPrognostic(IssmDouble kmax){
    336336
    337337        const int numdof=NUMVERTICES*NDOF1;
    338         double penalty_factor;
     338        IssmDouble penalty_factor;
    339339
    340340        /*Initialize Element vector and return if necessary*/
    … …  
    345345
    346346        //Create elementary matrix: add penalty to
    347         Ke->values[0*numdof+0]=+kmax*pow((double)10.0,penalty_factor);
    348         Ke->values[0*numdof+1]=-kmax*pow((double)10.0,penalty_factor);
    349         Ke->values[1*numdof+0]=-kmax*pow((double)10.0,penalty_factor);
    350         Ke->values[1*numdof+1]=+kmax*pow((double)10.0,penalty_factor);
     347        Ke->values[0*numdof+0]=+kmax*pow((IssmDouble)10.0,penalty_factor);
     348        Ke->values[0*numdof+1]=-kmax*pow((IssmDouble)10.0,penalty_factor);
     349        Ke->values[1*numdof+0]=-kmax*pow((IssmDouble)10.0,penalty_factor);
     350        Ke->values[1*numdof+1]=+kmax*pow((IssmDouble)10.0,penalty_factor);
    351351
    352352        /*Clean up and return*/

  • issm/trunk-jpl/src/c/objects/Loads/Penpair.h

    r12365 r12472  
    4141                /*}}}*/
    4242                /*Update virtual functions resolution: {{{*/
    43                 void  InputUpdateFromVector(double* vector, int name, int type);
     43                void  InputUpdateFromVector(IssmDouble* vector, int name, int type);
    4444                void  InputUpdateFromVector(int* vector, int name, int type);
    4545                void  InputUpdateFromVector(bool* vector, int name, int type);
    46                 void  InputUpdateFromMatrixDakota(double* matrix, int nrow, int ncols,int name, int type){_error_("Not implemented yet!");}
    47                 void  InputUpdateFromVectorDakota(double* vector, int name, int type){_error_("Not implemented yet!");}
     46                void  InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrow, int ncols,int name, int type){_error_("Not implemented yet!");}
     47                void  InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){_error_("Not implemented yet!");}
    4848                void  InputUpdateFromVectorDakota(int* vector, int name, int type){_error_("Not implemented yet!");}
    4949                void  InputUpdateFromVectorDakota(bool* vector, int name, int type){_error_("Not implemented yet!");}
    50                 void  InputUpdateFromConstant(double constant, int name);
     50                void  InputUpdateFromConstant(IssmDouble constant, int name);
    5151                void  InputUpdateFromConstant(int constant, int name);
    5252                void  InputUpdateFromConstant(bool constant, int name);
    53                 void  InputUpdateFromSolution(double* solution){_error_("Not implemented yet!");}
     53                void  InputUpdateFromSolution(IssmDouble* solution){_error_("Not implemented yet!");}
    5454                void  InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");};
    5555                /*}}}*/
    … …  
    6060                void  CreatePVector(Vector* pf);
    6161                void  CreateJacobianMatrix(Matrix* Jff);
    62                 void  PenaltyCreateKMatrix(Matrix* Kff,Matrix* Kfs,double kmax);
    63                 void  PenaltyCreatePVector(Vector* pf, double kmax);
    64                 void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax);
     62                void  PenaltyCreateKMatrix(Matrix* Kff,Matrix* Kfs,IssmDouble kmax);
     63                void  PenaltyCreatePVector(Vector* pf, IssmDouble kmax);
     64                void  PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax);
    6565                bool  InAnalysis(int analysis_type);
    6666                /*}}}*/
    6767                        /*Penpair management: {{{*/
    68                 ElementMatrix* PenaltyCreateKMatrixDiagnosticHoriz(double kmax);
    69                 ElementMatrix* PenaltyCreateKMatrixDiagnosticMacAyealPattyn(double kmax);
    70                 ElementMatrix* PenaltyCreateKMatrixDiagnosticStokes(double kmax);
    71                 ElementMatrix* PenaltyCreateKMatrixPrognostic(double kmax);
     68                ElementMatrix* PenaltyCreateKMatrixDiagnosticHoriz(IssmDouble kmax);
     69                ElementMatrix* PenaltyCreateKMatrixDiagnosticMacAyealPattyn(IssmDouble kmax);
     70                ElementMatrix* PenaltyCreateKMatrixDiagnosticStokes(IssmDouble kmax);
     71                ElementMatrix* PenaltyCreateKMatrixPrognostic(IssmDouble kmax);
    7272                /*}}}*/
    7373};

  • issm/trunk-jpl/src/c/objects/Loads/Riftfront.cpp

    r12365 r12472  
    4545        int    riftfront_type;
    4646        int    riftfront_fill;
    47         double riftfront_friction;
    48         double riftfront_fractionincrement;
     47        IssmDouble riftfront_friction;
     48        IssmDouble riftfront_fractionincrement;
    4949        bool   riftfront_shelf;
    5050        int    numberofelements;
    … …  
    135135        Input* input=NULL;
    136136        int fill;
    137         double friction,fractionincrement;
     137        IssmDouble friction,fractionincrement;
    138138
    139139       
    … …  
    258258}
    259259/*}}}*/
    260 /*FUNCTION Riftfront::InputUpdateFromConstant(double constant,int name) {{{*/
    261 void  Riftfront::InputUpdateFromConstant(double constant,int name){
     260/*FUNCTION Riftfront::InputUpdateFromConstant(IssmDouble constant,int name) {{{*/
     261void  Riftfront::InputUpdateFromConstant(IssmDouble constant,int name){
    262262
    263263        /*Check that name is a Riftfront input*/
    … …  
    269269}
    270270/*}}}*/
    271 /*FUNCTION Riftfront::InputUpdateFromConstant(double* constant,int name) {{{*/
    272 void    Riftfront::InputUpdateFromVector(double* vector, int name, int type){
     271/*FUNCTION Riftfront::InputUpdateFromConstant(IssmDouble* constant,int name) {{{*/
     272void    Riftfront::InputUpdateFromVector(IssmDouble* vector, int name, int type){
    273273
    274274        /*Check that name is a Riftfront input*/
    … …  
    309309/*}}}*/
    310310/*FUNCTION Riftfront::PenaltyCreateKMatrix {{{*/
    311 void  Riftfront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){
     311void  Riftfront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,IssmDouble kmax){
    312312
    313313        /*Retrieve parameters: */
    … …  
    335335/*}}}*/
    336336/*FUNCTION Riftfront::PenaltyCreatePVector {{{*/
    337 void  Riftfront::PenaltyCreatePVector(Vector* pf,double kmax){
     337void  Riftfront::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){
    338338
    339339        /*Retrieve parameters: */
    … …  
    381381/*Riftfront numerics*/
    382382/*FUNCTION Riftfront::PenaltyCreateKMatrixDiagnosticHoriz {{{*/
    383 ElementMatrix* Riftfront::PenaltyCreateKMatrixDiagnosticHoriz(double kmax){
     383ElementMatrix* Riftfront::PenaltyCreateKMatrixDiagnosticHoriz(IssmDouble kmax){
    384384
    385385        const int   numdof = NDOF2*NUMVERTICES;
    386386        int         i,j;
    387387        int         dofs[1]             = {0};
    388         double      Ke_gg[4][4];
    389         double      thickness;
    390         double      h[2];
    391         double      penalty_offset;
    392         double      friction;
     388        IssmDouble      Ke_gg[4][4];
     389        IssmDouble      thickness;
     390        IssmDouble      h[2];
     391        IssmDouble      penalty_offset;
     392        IssmDouble      friction;
    393393
    394394        /*Objects: */
    … …  
    464464/*}}}*/
    465465/*FUNCTION Riftfront::PenaltyCreatePVectorDiagnosticHoriz {{{*/
    466 ElementVector* Riftfront::PenaltyCreatePVectorDiagnosticHoriz(double kmax){
     466ElementVector* Riftfront::PenaltyCreatePVectorDiagnosticHoriz(IssmDouble kmax){
    467467
    468468        const int   numdof = NDOF2*NUMVERTICES;
    469469        int         i,j;
    470         double      rho_ice;
    471         double      rho_water;
    472         double      gravity;
    473         double      thickness;
    474         double      h[2];
    475         double      bed;
    476         double      b[2];
    477         double      pressure;
    478         double      pressure_litho;
    479         double      pressure_air;
    480         double      pressure_melange;
    481         double      pressure_water;
     470        IssmDouble      rho_ice;
     471        IssmDouble      rho_water;
     472        IssmDouble      gravity;
     473        IssmDouble      thickness;
     474        IssmDouble      h[2];
     475        IssmDouble      bed;
     476        IssmDouble      b[2];
     477        IssmDouble      pressure;
     478        IssmDouble      pressure_litho;
     479        IssmDouble      pressure_air;
     480        IssmDouble      pressure_melange;
     481        IssmDouble      pressure_water;
    482482        int         fill;
    483483        bool        shelf;
    … …  
    521521                if(shelf){
    522522                        /*We are on an ice shelf, hydrostatic equilibrium is used to determine the pressure for water fill: */
    523                         pressure=rho_ice*gravity*pow(thickness,(double)2)/(double)2  - rho_water*gravity*pow(bed,(double)2)/(double)2;
     523                        pressure=rho_ice*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2  - rho_water*gravity*pow(bed,(IssmDouble)2)/(IssmDouble)2;
    524524                }
    525525                else{
    526526                        //We are on an icesheet, we assume the water column fills the entire front: */
    527                         pressure=rho_ice*gravity*pow(thickness,(double)2)/(double)2  - rho_water*gravity*pow(thickness,(double)2)/(double)2;
     527                        pressure=rho_ice*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2  - rho_water*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2;
    528528                }
    529529        }
    530530        else if(fill==AirEnum){
    531                 pressure=rho_ice*gravity*pow(thickness,(double)2)/(double)2;   //icefront on an ice sheet, pressure imbalance ice vs air.
     531                pressure=rho_ice*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2;   //icefront on an ice sheet, pressure imbalance ice vs air.
    532532        }
    533533        else if(fill==IceEnum){ //icefront finding itself against another icefront (pressure imbalance is fully compensated, ice vs ice)
    … …  
    538538                if(!shelf) _error_("%s%i%s","fill type ",fill," not supported on ice sheets yet.");
    539539
    540                 pressure_litho=rho_ice*gravity*pow(thickness,(double)2)/(double)2;
     540                pressure_litho=rho_ice*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2;
    541541                pressure_air=0;
    542                 pressure_melange=rho_ice*gravity*pow(fraction*thickness,(double)2)/(double)2;
     542                pressure_melange=rho_ice*gravity*pow(fraction*thickness,(IssmDouble)2)/(IssmDouble)2;
    543543                pressure_water=1.0/2.0*rho_water*gravity*  ( pow(bed,2.0)-pow(rho_ice/rho_water*fraction*thickness,2.0) );
    544544
    … …  
    569569
    570570        const int   numnodes        = 2;
    571         double      max_penetration;
    572         double      penetration;
     571        IssmDouble      max_penetration;
     572        IssmDouble      penetration;
    573573        int         activate;
    574574        int         found;
    575575        int         unstable;
    576         double      vx1;
    577         double      vy1;
    578         double      vx2;
    579         double      vy2;
    580         double      fractionincrement;
     576        IssmDouble      vx1;
     577        IssmDouble      vy1;
     578        IssmDouble      vx2;
     579        IssmDouble      vy2;
     580        IssmDouble      fractionincrement;
    581581
    582582        /*Objects: */
    … …  
    632632                /*increase melange fraction: */
    633633                this->fraction+=fractionincrement;
    634                 if (this->fraction>1)this->fraction=(double)1.0;
     634                if (this->fraction>1)this->fraction=(IssmDouble)1.0;
    635635                //printf("riftfront %i fraction: %g\n",this->Id(),this->fraction);
    636636        }
    … …  
    675675
    676676        int found=0;
    677         double converged=0;
     677        IssmDouble converged=0;
    678678
    679679        this->inputs->GetInputValue(&converged,ConvergedEnum);
    … …  
    690690/*}}}*/
    691691/*FUNCTION Riftfront::MaxPenetration {{{*/
    692 int   Riftfront::MaxPenetration(double* ppenetration){
     692int   Riftfront::MaxPenetration(IssmDouble* ppenetration){
    693693
    694694        const int     numnodes=2;
    695         double        max_penetration;
    696         double        penetration=0;
     695        IssmDouble        max_penetration;
     696        IssmDouble        penetration=0;
    697697        int           found;
    698         double      vx1;
    699         double      vy1;
    700         double      vx2;
    701         double      vy2;
     698        IssmDouble      vx1;
     699        IssmDouble      vy1;
     700        IssmDouble      vx2;
     701        IssmDouble      vy2;
    702702
    703703        /*Objects: */
    … …  
    736736/*}}}*/
    737737/*FUNCTION Riftfront::Penetration {{{*/
    738 int   Riftfront::Penetration(double* ppenetration){
    739 
    740         double    vx1;
    741         double    vy1;
    742         double    vx2;
    743         double    vy2;
    744 
    745         double    penetration;
     738int   Riftfront::Penetration(IssmDouble* ppenetration){
     739
     740        IssmDouble    vx1;
     741        IssmDouble    vy1;
     742        IssmDouble    vx2;
     743        IssmDouble    vy2;
     744
     745        IssmDouble    penetration;
    746746        int       found;
    747747
    … …  
    779779
    780780        const int   numnodes        = 2;
    781         double      max_penetration;
    782         double      penetration;
     781        IssmDouble      max_penetration;
     782        IssmDouble      penetration;
    783783        int         activate;
    784784        int         unstable;
    785785        int         found;
    786         double      vx1;
    787         double      vy1;
    788         double      vx2;
    789         double      vy2;
     786        IssmDouble      vx1;
     787        IssmDouble      vy1;
     788        IssmDouble      vx2;
     789        IssmDouble      vy2;
    790790
    791791        /*Objects: */
    … …  
    831831
    832832        const int   numnodes    = 2;
    833         double      penetration;
     833        IssmDouble      penetration;
    834834        int         unstable;
    835835        int         found;
    836         double      vx1;
    837         double      vy1;
    838         double      vx2;
    839         double      vy2;
     836        IssmDouble      vx1;
     837        IssmDouble      vy1;
     838        IssmDouble      vx2;
     839        IssmDouble      vy2;
    840840
    841841        /*Objects: */

  • issm/trunk-jpl/src/c/objects/Loads/Riftfront.h

    r12365 r12472  
    3838                bool     prestable;
    3939                bool     material_converged;
    40                 double   normal[2];
    41                 double   length;
    42                 double   fraction;
     40                IssmDouble   normal[2];
     41                IssmDouble   length;
     42                IssmDouble   fraction;
    4343                int      state;
    4444               
    … …  
    6161                /*}}}*/
    6262                /*Update virtual functions resolution: {{{*/
    63                 void    InputUpdateFromVector(double* vector, int name, int type);
     63                void    InputUpdateFromVector(IssmDouble* vector, int name, int type);
    6464                void    InputUpdateFromVector(int* vector, int name, int type){_error_("Not implemented yet!");}
    6565                void    InputUpdateFromVector(bool* vector, int name, int type){_error_("Not implemented yet!");}
    66                 void    InputUpdateFromMatrixDakota(double* matrix, int nrows,int ncols, int name, int type){_error_("Not implemented yet!");}
    67                 void    InputUpdateFromVectorDakota(double* vector, int name, int type){_error_("Not implemented yet!");}
     66                void    InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows,int ncols, int name, int type){_error_("Not implemented yet!");}
     67                void    InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){_error_("Not implemented yet!");}
    6868                void    InputUpdateFromVectorDakota(int* vector, int name, int type){_error_("Not implemented yet!");}
    6969                void    InputUpdateFromVectorDakota(bool* vector, int name, int type){_error_("Not implemented yet!");}
    70                 void    InputUpdateFromConstant(double constant, int name);
     70                void    InputUpdateFromConstant(IssmDouble constant, int name);
    7171                void    InputUpdateFromConstant(int constant, int name){_error_("Not implemented yet!");}
    7272                void    InputUpdateFromConstant(bool constant, int name);
    73                 void    InputUpdateFromSolution(double* solution){_error_("Not implemented yet!");}
     73                void    InputUpdateFromSolution(IssmDouble* solution){_error_("Not implemented yet!");}
    7474                void  InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");};
    7575                /*}}}*/
    … …  
    8080                void  CreatePVector(Vector* pf);
    8181                void  CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");};
    82                 void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");};
    83                 void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax);
    84                 void  PenaltyCreatePVector(Vector* pf, double kmax);
     82                void  PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){_error_("Not implemented yet");};
     83                void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, IssmDouble kmax);
     84                void  PenaltyCreatePVector(Vector* pf, IssmDouble kmax);
    8585                bool  InAnalysis(int analysis_type);
    8686                /*}}}*/
    8787                /*Riftfront specific routines: {{{*/
    8888                bool  PreStable();
    89                 ElementMatrix* PenaltyCreateKMatrixDiagnosticHoriz(double kmax);
    90                 ElementVector* PenaltyCreatePVectorDiagnosticHoriz(double kmax);
     89                ElementMatrix* PenaltyCreateKMatrixDiagnosticHoriz(IssmDouble kmax);
     90                ElementVector* PenaltyCreatePVectorDiagnosticHoriz(IssmDouble kmax);
    9191                void  SetPreStable();
    9292                int   PreConstrain(int* punstable);
    … …  
    9494                void  FreezeConstraints(void);
    9595                bool  IsFrozen(void);
    96                 int   Penetration(double* ppenetration);
    97                 int   MaxPenetration(double* ppenetration);
     96                int   Penetration(IssmDouble* ppenetration);
     97                int   MaxPenetration(IssmDouble* ppenetration);
    9898                int   PotentialUnstableConstraint(int* punstable);
    9999                int   IsMaterialStable(void);

  • issm/trunk-jpl/src/c/objects/Materials/Matice.cpp

    r12365 r12472  
    129129/*}}}*/
    130130/*FUNCTION Matice::GetB {{{*/
    131 double Matice::GetB(){
     131IssmDouble Matice::GetB(){
    132132
    133133        /*Output*/
    134         double B;
     134        IssmDouble B;
    135135
    136136        inputs->GetInputAverage(&B,MaterialsRheologyBEnum);
    … …  
    139139/*}}}*/
    140140/*FUNCTION Matice::GetBbar {{{*/
    141 double Matice::GetBbar(){
     141IssmDouble Matice::GetBbar(){
    142142
    143143        /*Output*/
    144         double Bbar;
     144        IssmDouble Bbar;
    145145
    146146        inputs->GetInputAverage(&Bbar,MaterialsRheologyBbarEnum);
    … …  
    149149/*}}}*/
    150150/*FUNCTION Matice::GetN {{{*/
    151 double Matice::GetN(){
     151IssmDouble Matice::GetN(){
    152152
    153153        /*Output*/
    154         double n;
     154        IssmDouble n;
    155155
    156156        inputs->GetInputAverage(&n,MaterialsRheologyNEnum);
    … …  
    191191/*}}}*/
    192192/*FUNCTION Matice::GetViscosity2d {{{*/
    193 void  Matice::GetViscosity2d(double* pviscosity, double* epsilon){
     193void  Matice::GetViscosity2d(IssmDouble* pviscosity, IssmDouble* epsilon){
    194194        /*From a string tensor and a material object, return viscosity, using Glen's flow law.
    195195                                                                                                    B
    … …  
    205205
    206206        /*output: */
    207         double viscosity;
     207        IssmDouble viscosity;
    208208
    209209        /*input strain rate: */
    210         double exx,eyy,exy;
     210        IssmDouble exx,eyy,exy;
    211211
    212212        /*Intermediary: */
    213         double A,e;
    214         double B,n;
     213        IssmDouble A,e;
     214        IssmDouble B,n;
    215215
    216216        /*Get B and n*/
    … …  
    224224        else{
    225225                if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0)){
    226                         viscosity=0.5*pow((double)10,(double)14);
     226                        viscosity=0.5*pow((IssmDouble)10,(IssmDouble)14);
    227227                }
    228228                else{
    … …  
    236236                        if(A==0){
    237237                                /*Maxiviscositym viscosity for 0 shear areas: */
    238                                 viscosity=2.5*pow((double)10,(double)17);
     238                                viscosity=2.5*pow((IssmDouble)10,(IssmDouble)17);
    239239                        }
    240240                        else{
    … …  
    255255/*}}}*/
    256256/*FUNCTION Matice::GetViscosity3d {{{*/
    257 void  Matice::GetViscosity3d(double* pviscosity3d, double* epsilon){
     257void  Matice::GetViscosity3d(IssmDouble* pviscosity3d, IssmDouble* epsilon){
    258258
    259259        /*Return viscosity accounting for steady state power law creep [Thomas and MacAyeal, 1982]:
    … …  
    271271       
    272272        /*output: */
    273         double viscosity3d;
     273        IssmDouble viscosity3d;
    274274
    275275        /*input strain rate: */
    276         double exx,eyy,exy,exz,eyz;
     276        IssmDouble exx,eyy,exy,exz,eyz;
    277277
    278278        /*Intermediaries: */
    279         double A,e;
    280         double B,n;
     279        IssmDouble A,e;
     280        IssmDouble B,n;
    281281
    282282        /*Get B and n*/
    … …  
    291291                if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0) &&
    292292                                (epsilon[3]==0) && (epsilon[4]==0)){
    293                         viscosity3d=0.5*pow((double)10,(double)14);
     293                        viscosity3d=0.5*pow((IssmDouble)10,(IssmDouble)14);
    294294                }
    295295                else{
    … …  
    306306                        if(A==0){
    307307                                /*Maxiviscosity3dm viscosity for 0 shear areas: */
    308                                 viscosity3d=2.25*pow((double)10,(double)17);
     308                                viscosity3d=2.25*pow((IssmDouble)10,(IssmDouble)17);
    309309                        }
    310310                        else{
    … …  
    326326/*}}}*/
    327327/*FUNCTION Matice::GetViscosity3dStokes {{{*/
    328 void  Matice::GetViscosity3dStokes(double* pviscosity3d, double* epsilon){
     328void  Matice::GetViscosity3dStokes(IssmDouble* pviscosity3d, IssmDouble* epsilon){
    329329        /*Return viscosity accounting for steady state power law creep [Thomas and MacAyeal, 1982]:
    330330         *
    … …  
    341341       
    342342        /*output: */
    343         double viscosity3d;
     343        IssmDouble viscosity3d;
    344344
    345345        /*input strain rate: */
    346         double exx,eyy,exy,exz,eyz,ezz;
     346        IssmDouble exx,eyy,exy,exz,eyz,ezz;
    347347
    348348        /*Intermediaries: */
    349         double A,e;
    350         double B,n;
    351         double eps0;
     349        IssmDouble A,e;
     350        IssmDouble B,n;
     351        IssmDouble eps0;
    352352
    353353        /*Get B and n*/
    354         eps0=pow((double)10,(double)-27);
     354        eps0=pow((IssmDouble)10,(IssmDouble)-27);
    355355        B=GetB();
    356356        n=GetN();
    … …  
    363363                if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0) &&
    364364                                (epsilon[3]==0) && (epsilon[4]==0) && (epsilon[5]==0)){
    365                         viscosity3d=0.5*pow((double)10,(double)14);
     365                        viscosity3d=0.5*pow((IssmDouble)10,(IssmDouble)14);
    366366                }
    367367                else{
    … …  
    379379                        if(A==0){
    380380                                /*Maxiviscosity3dm viscosity for 0 shear areas: */
    381                                 viscosity3d=2.25*pow((double)10,(double)17);
     381                                viscosity3d=2.25*pow((IssmDouble)10,(IssmDouble)17);
    382382                        }
    383383                        else{
    … …  
    398398/*}}}*/
    399399/*FUNCTION Matice::GetViscosityComplement {{{*/
    400 void  Matice::GetViscosityComplement(double* pviscosity_complement, double* epsilon){
     400void  Matice::GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* epsilon){
    401401        /*Return viscosity accounting for steady state power law creep [Thomas and MacAyeal, 1982]:
    402402         *
    … …  
    410410       
    411411        /*output: */
    412         double viscosity_complement;
     412        IssmDouble viscosity_complement;
    413413
    414414        /*input strain rate: */
    415         double exx,eyy,exy;
     415        IssmDouble exx,eyy,exy;
    416416
    417417        /*Intermediary value A and exponent e: */
    418         double A,e;
    419         double B,n;
     418        IssmDouble A,e;
     419        IssmDouble B,n;
    420420
    421421        /*Get B and n*/
    … …  
    432432                if(A==0){
    433433                        /*Maximum viscosity_complement for 0 shear areas: */
    434                         viscosity_complement=2.25*pow((double)10,(double)17);
     434                        viscosity_complement=2.25*pow((IssmDouble)10,(IssmDouble)17);
    435435                }
    436436                else{
    … …  
    441441        }
    442442        else{
    443                 viscosity_complement=4.5*pow((double)10,(double)17);
     443                viscosity_complement=4.5*pow((IssmDouble)10,(IssmDouble)17);
    444444        }
    445445
    … …  
    454454/*}}}*/
    455455/*FUNCTION Matice::GetViscosityDerivativeEpsSquare{{{*/
    456 void  Matice::GetViscosityDerivativeEpsSquare(double* pmu_prime, double* epsilon){
     456void  Matice::GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* epsilon){
    457457
    458458        /*output: */
    459         double mu_prime;
    460         double mu,n,eff2;
     459        IssmDouble mu_prime;
     460        IssmDouble mu,n,eff2;
    461461
    462462        /*input strain rate: */
    463         double exx,eyy,exy,exz,eyz;
     463        IssmDouble exx,eyy,exy,exz,eyz;
    464464
    465465        /*Get visocisty and n*/
    … …  
    469469        if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0) &&
    470470                                (epsilon[3]==0) && (epsilon[4]==0)){
    471                 mu_prime=0.5*pow((double)10,(double)14);
     471                mu_prime=0.5*pow((IssmDouble)10,(IssmDouble)14);
    472472        }
    473473        else{
    … …  
    488488/*}}}*/
    489489/*FUNCTION Matice::GetViscosity2dDerivativeEpsSquare{{{*/
    490 void  Matice::GetViscosity2dDerivativeEpsSquare(double* pmu_prime, double* epsilon){
     490void  Matice::GetViscosity2dDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* epsilon){
    491491
    492492        /*output: */
    493         double mu_prime;
    494         double mu,n,eff2;
     493        IssmDouble mu_prime;
     494        IssmDouble mu,n,eff2;
    495495
    496496        /*input strain rate: */
    497         double exx,eyy,exy,exz;
     497        IssmDouble exx,eyy,exy,exz;
    498498
    499499        /*Get visocisty and n*/
    … …  
    502502
    503503        if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0)){
    504                 mu_prime=0.5*pow((double)10,(double)14);
     504                mu_prime=0.5*pow((IssmDouble)10,(IssmDouble)14);
    505505        }
    506506        else{
    … …  
    526526}
    527527/*}}}*/
    528 /*FUNCTION Matice::InputUpdateFromVector(double* vector, int name, int type) {{{*/
    529 void  Matice::InputUpdateFromVector(double* vector, int name, int type){
     528/*FUNCTION Matice::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/
     529void  Matice::InputUpdateFromVector(IssmDouble* vector, int name, int type){
    530530
    531531        /*Intermediaries*/
    … …  
    545545
    546546                                case TriaEnum:
    547                                         double values[3];
     547                                        IssmDouble values[3];
    548548                                        for (int i=0;i<3;i++) values[i]=vector[((Tria*)element)->nodes[i]->GetVertexDof()];
    549549                                        this->inputs->AddInput(new TriaP1Input(name,values));
    … …  
    566566}
    567567/*}}}*/
    568 /*FUNCTION Matice::InputUpdateFromVectorDakota(double* vector, int name, int type) {{{*/
    569 void  Matice::InputUpdateFromVectorDakota(double* vector, int name, int type){
     568/*FUNCTION Matice::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type) {{{*/
     569void  Matice::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){
    570570
    571571        /*Intermediaries*/
    … …  
    587587
    588588                                case TriaEnum:
    589                                         double values[3];
     589                                        IssmDouble values[3];
    590590                                        for (int i=0;i<3;i++) values[i]=vector[((Tria*)element)->nodes[i]->GetSidList()]; //use sid list, to index into serial oriented vector
    591591                                        this->inputs->AddInput(new TriaP1Input(name,values));
    … …  
    618618/*}}}*/
    619619/*FUNCTION Matice::InputUpdateFromMatrixDakota(int* vector, int name, int type) {{{*/
    620 void  Matice::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols,int name, int type){
     620void  Matice::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols,int name, int type){
    621621        /*Nothing updated yet*/
    622622}
    … …  
    632632}
    633633/*}}}*/
    634 /*FUNCTION Matice::InputUpdateFromConstant(double constant, int name) {{{*/
    635 void  Matice::InputUpdateFromConstant(double constant, int name){
     634/*FUNCTION Matice::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/
     635void  Matice::InputUpdateFromConstant(IssmDouble constant, int name){
    636636        /*Nothing updated yet*/
    637637}
    … …  
    648648/*}}}*/
    649649/*FUNCTION Matice::InputUpdateFromSolution{{{*/
    650 void  Matice::InputUpdateFromSolution(double* solution){
     650void  Matice::InputUpdateFromSolution(IssmDouble* solution){
    651651        /*Nothing updated yet*/
    652652}
    … …  
    671671                /*Intermediaries*/
    672672                const int num_vertices = 3; //Tria has 3 vertices
    673                 double    nodeinputs[num_vertices];
    674                 double    cmmininputs[num_vertices];
    675                 double    cmmaxinputs[num_vertices];
     673                IssmDouble    nodeinputs[num_vertices];
     674                IssmDouble    cmmininputs[num_vertices];
     675                IssmDouble    cmmaxinputs[num_vertices];
    676676
    677677                /*Get B*/
    … …  
    713713                /*Intermediaries*/
    714714                const int num_vertices = 6; //Penta has 6 vertices
    715                 double    nodeinputs[num_vertices];
    716                 double    cmmininputs[num_vertices];
    717                 double    cmmaxinputs[num_vertices];
     715                IssmDouble    nodeinputs[num_vertices];
     716                IssmDouble    cmmininputs[num_vertices];
     717                IssmDouble    cmmaxinputs[num_vertices];
    718718
    719719                /*Get B*/

  • issm/trunk-jpl/src/c/objects/Materials/Matice.h

    r12365 r12472  
    3939                /*}}}*/
    4040                /*Update virtual functions definitions: {{{*/
    41                 void  InputUpdateFromVector(double* vector, int name, int type);
     41                void  InputUpdateFromVector(IssmDouble* vector, int name, int type);
    4242                void  InputUpdateFromVector(int* vector, int name, int type);
    4343                void  InputUpdateFromVector(bool* vector, int name, int type);
    44                 void  InputUpdateFromMatrixDakota(double* matrix, int nrow, int ncols, int name, int type);
    45                 void  InputUpdateFromVectorDakota(double* vector, int name, int type);
     44                void  InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrow, int ncols, int name, int type);
     45                void  InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type);
    4646                void  InputUpdateFromVectorDakota(int* vector, int name, int type);
    4747                void  InputUpdateFromVectorDakota(bool* vector, int name, int type);
    48                 void  InputUpdateFromConstant(double constant, int name);
     48                void  InputUpdateFromConstant(IssmDouble constant, int name);
    4949                void  InputUpdateFromConstant(int constant, int name);
    5050                void  InputUpdateFromConstant(bool constant, int name);
    51                 void  InputUpdateFromSolution(double* solution);
     51                void  InputUpdateFromSolution(IssmDouble* solution);
    5252                void  InputUpdateFromIoModel(int index, IoModel* iomodel);
    5353                /*}}}*/
    … …  
    5959                /*Matice Numerics: {{{*/
    6060                void   SetCurrentConfiguration(Elements* elementsin,Loads* loadsin,Nodes* nodesin,Vertices* verticesin,Materials* materialsin,Parameters* parametersin);
    61                 void   GetViscosity2d(double* pviscosity, double* pepsilon);
    62                 void   GetViscosity3d(double* pviscosity3d, double* pepsilon);
    63                 void   GetViscosity3dStokes(double* pviscosity3d, double* epsilon);
    64                 void   GetViscosityComplement(double* pviscosity_complement, double* pepsilon);
    65                 void   GetViscosityDerivativeEpsSquare(double* pmu_prime, double* pepsilon);
    66                 void   GetViscosity2dDerivativeEpsSquare(double* pmu_prime, double* pepsilon);
    67                 double GetB();
    68                 double GetBbar();
    69                 double GetN();
     61                void   GetViscosity2d(IssmDouble* pviscosity, IssmDouble* pepsilon);
     62                void   GetViscosity3d(IssmDouble* pviscosity3d, IssmDouble* pepsilon);
     63                void   GetViscosity3dStokes(IssmDouble* pviscosity3d, IssmDouble* epsilon);
     64                void   GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
     65                void   GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
     66                void   GetViscosity2dDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
     67                IssmDouble GetB();
     68                IssmDouble GetBbar();
     69                IssmDouble GetN();
    7070                bool   IsInput(int name);
    7171                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Materials/Matpar.cpp

    r12365 r12472  
    104104
    105105/*Update virtual functions definitions:*/
    106 /*FUNCTION Matpar::InputUpdateFromVector(double* vector, int name, int type) {{{*/
    107 void   Matpar::InputUpdateFromVector(double* vector, int name, int type){
     106/*FUNCTION Matpar::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/
     107void   Matpar::InputUpdateFromVector(IssmDouble* vector, int name, int type){
    108108        /*Nothing updated yet*/
    109109}
    … …  
    119119}
    120120/*}}}*/
    121 /*FUNCTION Matpar::InputUpdateFromVectorDakota(double* vector, int name, int type) {{{*/
    122 void   Matpar::InputUpdateFromVectorDakota(double* vector, int name, int type){
     121/*FUNCTION Matpar::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type) {{{*/
     122void   Matpar::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){
    123123        /*Nothing updated yet*/
    124124}
    … …  
    135135/*}}}*/
    136136/*FUNCTION Matpar::InputUpdateFromMatrixDakota(int* vector, int name, int type) {{{*/
    137 void  Matpar::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols,int name, int type){
    138         /*Nothing updated yet*/
    139 }
    140 /*}}}*/
    141 /*FUNCTION Matpar::InputUpdateFromConstant(double constant, int name) {{{*/
    142 void   Matpar::InputUpdateFromConstant(double constant, int name){
     137void  Matpar::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols,int name, int type){
     138        /*Nothing updated yet*/
     139}
     140/*}}}*/
     141/*FUNCTION Matpar::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/
     142void   Matpar::InputUpdateFromConstant(IssmDouble constant, int name){
    143143
    144144        switch(name){
    … …  
    199199/*}}}*/
    200200/*FUNCTION Matpar::InputUpdateFromSolution{{{*/
    201 void   Matpar::InputUpdateFromSolution(double* solution){
     201void   Matpar::InputUpdateFromSolution(IssmDouble* solution){
    202202        /*Nothing updated yet*/
    203203}
    … …  
    213213/*}}}*/
    214214/*FUNCTION Matpar::GetBeta {{{*/
    215 double Matpar::GetBeta(){
     215IssmDouble Matpar::GetBeta(){
    216216        return beta;
    217217}
    218218/*}}}*/
    219219/*FUNCTION Matpar::GetG {{{*/
    220 double Matpar::GetG(){
     220IssmDouble Matpar::GetG(){
    221221        return g;
    222222}
    223223/*}}}*/
    224224/*FUNCTION Matpar::GetHeatCapacity {{{*/
    225 double Matpar::GetHeatCapacity(){
     225IssmDouble Matpar::GetHeatCapacity(){
    226226        return heatcapacity;
    227227}
    228228/*}}}*/
    229229/*FUNCTION Matpar::GetLatentHeat {{{*/
    230 double Matpar::GetLatentHeat(){
     230IssmDouble Matpar::GetLatentHeat(){
    231231        return latentheat;
    232232}
    233233/*}}}*/
    234234/*FUNCTION Matpar::GetMeltingPoint {{{*/
    235 double Matpar::GetMeltingPoint(){
     235IssmDouble Matpar::GetMeltingPoint(){
    236236        return meltingpoint;
    237237}
    238238/*}}}*/
    239239/*FUNCTION Matpar::GetReferenceTemperature {{{*/
    240 double Matpar::GetReferenceTemperature(){
     240IssmDouble Matpar::GetReferenceTemperature(){
    241241        return referencetemperature;
    242242}
    243243/*}}}*/
    244244/*FUNCTION Matpar::GetMixedLayerCapacity {{{*/
    245 double Matpar::GetMixedLayerCapacity(){
     245IssmDouble Matpar::GetMixedLayerCapacity(){
    246246        return mixed_layer_capacity;
    247247}
    248248/*}}}*/
    249249/*FUNCTION Matpar::GetRhoIce {{{*/
    250 double Matpar::GetRhoIce(){
     250IssmDouble Matpar::GetRhoIce(){
    251251       
    252252        return rho_ice;
    … …  
    254254/*}}}*/
    255255/*FUNCTION Matpar::GetRhoWater {{{*/
    256 double Matpar::GetRhoWater(){
     256IssmDouble Matpar::GetRhoWater(){
    257257        return rho_water;
    258258}
    259259/*}}}*/
    260260/*FUNCTION Matpar::GetRhoFreshwater {{{*/
    261 double Matpar::GetRhoFreshwater(){
     261IssmDouble Matpar::GetRhoFreshwater(){
    262262        return rho_freshwater;
    263263}
    264264/*}}}*/
    265265/*FUNCTION Matpar::GetMuWater {{{*/
    266 double Matpar::GetMuWater(){
     266IssmDouble Matpar::GetMuWater(){
    267267        return mu_water;
    268268}
    269269/*}}}*/
    270270/*FUNCTION Matpar::GetThermalConductivity {{{*/
    271 double Matpar::GetThermalConductivity(){
     271IssmDouble Matpar::GetThermalConductivity(){
    272272        return thermalconductivity;
    273273}
    274274/*}}}*/
    275275/*FUNCTION Matpar::GetThermalExchangeVelocity {{{*/
    276 double Matpar::GetThermalExchangeVelocity(){
     276IssmDouble Matpar::GetThermalExchangeVelocity(){
    277277        return thermal_exchange_velocity;
    278278}
    279279/*}}}*/
    280280/*FUNCTION Matpar::GetKn {{{*/           
    281 double Matpar::GetKn(){                 
     281IssmDouble Matpar::GetKn(){                     
    282282        return kn;               
    283283}               
    284284/*}}}*/                 
    285285/*FUNCTION Matpar::GetHydrologyP {{{*/                   
    286 double Matpar::GetHydrologyP(){         
     286IssmDouble Matpar::GetHydrologyP(){             
    287287        return hydro_p;                 
    288288}               
    289289/*}}}*/                 
    290290/*FUNCTION Matqar::GetHydrologyQ {{{*/                   
    291 double Matpar::GetHydrologyQ(){         
     291IssmDouble Matpar::GetHydrologyQ(){             
    292292        return hydro_q;                 
    293293}               
    294294/*}}}*/                 
    295295/*FUNCTION Matpar::GetHydrologyCR {{{*/         
    296 double Matpar::GetHydrologyCR(){                 
     296IssmDouble Matpar::GetHydrologyCR(){             
    297297        return hydro_CR;                 
    298298}               
    299299/*}}}*/                 
    300300/*FUNCTION Matpar::GetHydrologyN {{{*/                   
    301 double Matpar::GetHydrologyN(){         
     301IssmDouble Matpar::GetHydrologyN(){             
    302302        return hydro_n;                 
    303303}               
    304304/*}}}*/
    305305/*FUNCTION Matpar::TMeltingPoint {{{*/
    306 double Matpar::TMeltingPoint(double pressure){
     306IssmDouble Matpar::TMeltingPoint(IssmDouble pressure){
    307307        return meltingpoint-beta*pressure;
    308308}
    309309/*}}}*/
    310310/*FUNCTION Matpar::PureIceEnthalpy{{{*/
    311 double Matpar::PureIceEnthalpy(double pressure){
     311IssmDouble Matpar::PureIceEnthalpy(IssmDouble pressure){
    312312        return heatcapacity*(TMeltingPoint(pressure)-referencetemperature);
    313313}
    314314/*}}}*/
    315315/*FUNCTION Matpar::GetEnthalpyDiffusionParameter{{{*/
    316 double Matpar::GetEnthalpyDiffusionParameter(double enthalpy,double pressure){
     316IssmDouble Matpar::GetEnthalpyDiffusionParameter(IssmDouble enthalpy,IssmDouble pressure){
    317317        if(enthalpy<PureIceEnthalpy(pressure)){
    318318                return thermalconductivity/(rho_ice*heatcapacity);
    … …  
    324324/*}}}*/
    325325/*FUNCTION Matpar::EnthalpyToThermal {{{*/
    326 void Matpar::EnthalpyToThermal(double* ptemperature,double* pwaterfraction,double enthalpy,double pressure){
     326void Matpar::EnthalpyToThermal(IssmDouble* ptemperature,IssmDouble* pwaterfraction,IssmDouble enthalpy,IssmDouble pressure){
    327327
    328328        /*Ouput*/
    329         double temperature,waterfraction;
     329        IssmDouble temperature,waterfraction;
    330330       
    331331        if(enthalpy<PureIceEnthalpy(pressure)){
    … …  
    344344/*}}}*/
    345345/*FUNCTION Matpar::ThermalToEnthalpy {{{*/
    346 void Matpar::ThermalToEnthalpy(double * penthalpy,double temperature,double waterfraction,double pressure){
     346void Matpar::ThermalToEnthalpy(IssmDouble * penthalpy,IssmDouble temperature,IssmDouble waterfraction,IssmDouble pressure){
    347347
    348348        /*Ouput*/
    349         double enthalpy;
     349        IssmDouble enthalpy;
    350350       
    351351        if(temperature<TMeltingPoint(pressure)){

  • issm/trunk-jpl/src/c/objects/Materials/Matpar.h

    r12365 r12472  
    1616        private:
    1717                int       mid;
    18                 double  rho_ice;
    19                 double  rho_water;
    20                 double  rho_freshwater;
    21                 double  mu_water;
    22                 double  heatcapacity;
    23                 double  thermalconductivity;
    24                 double  latentheat;
    25                 double  beta;
    26                 double  meltingpoint;
    27                 double  referencetemperature;
    28                 double  mixed_layer_capacity;
    29                 double  thermal_exchange_velocity;
    30                 double  g;
     18                IssmDouble  rho_ice;
     19                IssmDouble  rho_water;
     20                IssmDouble  rho_freshwater;
     21                IssmDouble  mu_water;
     22                IssmDouble  heatcapacity;
     23                IssmDouble  thermalconductivity;
     24                IssmDouble  latentheat;
     25                IssmDouble  beta;
     26                IssmDouble  meltingpoint;
     27                IssmDouble  referencetemperature;
     28                IssmDouble  mixed_layer_capacity;
     29                IssmDouble  thermal_exchange_velocity;
     30                IssmDouble  g;
    3131
    3232                /*hydrology: */         
    33                 double  kn;                     
    34                 double  hydro_p;                 
    35                 double  hydro_q;                 
    36                 double  hydro_CR;                       
    37                 double  hydro_n;
     33                IssmDouble  kn;                 
     34                IssmDouble  hydro_p;             
     35                IssmDouble  hydro_q;             
     36                IssmDouble  hydro_CR;                   
     37                IssmDouble  hydro_n;
    3838
    3939        public:
    … …  
    5151                /*}}}*/
    5252                /*Update virtual functions resolution: {{{*/
    53                 void   InputUpdateFromVector(double* vector, int name, int type);
     53                void   InputUpdateFromVector(IssmDouble* vector, int name, int type);
    5454                void   InputUpdateFromVector(int* vector, int name, int type);
    5555                void   InputUpdateFromVector(bool* vector, int name, int type);
    56                 void   InputUpdateFromMatrixDakota(double* matrix,int nrows,int ncols, int name, int type);
    57                 void   InputUpdateFromVectorDakota(double* vector, int name, int type);
     56                void   InputUpdateFromMatrixDakota(IssmDouble* matrix,int nrows,int ncols, int name, int type);
     57                void   InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type);
    5858                void   InputUpdateFromVectorDakota(int* vector, int name, int type);
    5959                void   InputUpdateFromVectorDakota(bool* vector, int name, int type);
    60                 void   InputUpdateFromConstant(double constant, int name);
     60                void   InputUpdateFromConstant(IssmDouble constant, int name);
    6161                void   InputUpdateFromConstant(int constant, int name);
    6262                void   InputUpdateFromConstant(bool constant, int name);
    63                 void   InputUpdateFromSolution(double* solution);
     63                void   InputUpdateFromSolution(IssmDouble* solution);
    6464                void   InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");};
    6565                /*}}}*/
    … …  
    7070                /*}}}*/
    7171                /*Numerics: {{{*/
    72                 double GetG();
    73                 double GetRhoIce();
    74                 double GetRhoWater();
    75                 double GetRhoFreshwater();
    76                 double GetMuWater();
    77                 double GetMixedLayerCapacity();
    78                 double GetThermalExchangeVelocity();
    79                 double GetHeatCapacity();
    80                 double GetThermalConductivity();
    81                 double GetLatentHeat();
    82                 double GetBeta();
    83                 double GetMeltingPoint();
    84                 double GetReferenceTemperature();
    85                 double GetKn();
    86                 double GetHydrologyP();
    87                 double GetHydrologyQ();
    88                 double GetHydrologyCR();
    89                 double GetHydrologyN();
    90                 double TMeltingPoint(double pressure);
    91                 double PureIceEnthalpy(double pressure);
    92                 double GetEnthalpyDiffusionParameter(double enthalpy,double pressure);
    93                 void   EnthalpyToThermal(double* ptemperature,double* pwaterfraction,double enthalpy,double pressure);
    94                 void   ThermalToEnthalpy(double* penthalpy,double temperature,double waterfraction,double pressure);
     72                IssmDouble GetG();
     73                IssmDouble GetRhoIce();
     74                IssmDouble GetRhoWater();
     75                IssmDouble GetRhoFreshwater();
     76                IssmDouble GetMuWater();
     77                IssmDouble GetMixedLayerCapacity();
     78                IssmDouble GetThermalExchangeVelocity();
     79                IssmDouble GetHeatCapacity();
     80                IssmDouble GetThermalConductivity();
     81                IssmDouble GetLatentHeat();
     82                IssmDouble GetBeta();
     83                IssmDouble GetMeltingPoint();
     84                IssmDouble GetReferenceTemperature();
     85                IssmDouble GetKn();
     86                IssmDouble GetHydrologyP();
     87                IssmDouble GetHydrologyQ();
     88                IssmDouble GetHydrologyCR();
     89                IssmDouble GetHydrologyN();
     90                IssmDouble TMeltingPoint(IssmDouble pressure);
     91                IssmDouble PureIceEnthalpy(IssmDouble pressure);
     92                IssmDouble GetEnthalpyDiffusionParameter(IssmDouble enthalpy,IssmDouble pressure);
     93                void   EnthalpyToThermal(IssmDouble* ptemperature,IssmDouble* pwaterfraction,IssmDouble enthalpy,IssmDouble pressure);
     94                void   ThermalToEnthalpy(IssmDouble* penthalpy,IssmDouble temperature,IssmDouble waterfraction,IssmDouble pressure);
    9595                /*}}}*/
    9696

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

    r12446 r12472  
    101101        /*Gset and values*/
    102102        this->gglobaldoflist=xNew<int>(this->nrows);
    103         this->values=xNewZeroInit<double>(this->nrows*this->ncols);
     103        this->values=xNewZeroInit<IssmDouble>(this->nrows*this->ncols);
    104104        for(i=0;i<Ke1->nrows;i++){
    105105                for(j=0;j<Ke1->ncols;j++){
    … …  
    205205
    206206        /*fill values with 0: */
    207         this->values=xNewZeroInit<double>(this->nrows*this->ncols);
     207        this->values=xNewZeroInit<IssmDouble>(this->nrows*this->ncols);
    208208
    209209        /*g list*/
    … …  
    231231ElementMatrix::~ElementMatrix(){
    232232       
    233         xDelete<double>(this->values);
     233        xDelete<IssmDouble>(this->values);
    234234        xDelete<int>(this->gglobaldoflist);
    235235        xDelete<int>(this->row_flocaldoflist);
    … …  
    249249
    250250        int i,j;
    251         double* localvalues=NULL;
     251        IssmDouble* localvalues=NULL;
    252252
    253253        /*If Kfs is not provided, call the other function*/
    … …  
    265265                if(this->row_fsize){
    266266                        /*first, retrieve values that are in the f-set from the g-set values matrix: */
    267                         localvalues=xNew<double>(this->row_fsize*this->row_fsize);
     267                        localvalues=xNew<IssmDouble>(this->row_fsize*this->row_fsize);
    268268                        for(i=0;i<this->row_fsize;i++){
    269269                                for(j=0;j<this->row_fsize;j++){
    … …  
    275275
    276276                        /*Free ressources:*/
    277                         xDelete<double>(localvalues);
     277                        xDelete<IssmDouble>(localvalues);
    278278                }
    279279
    … …  
    281281                if((this->row_ssize!=0) && (this->row_fsize!=0)){
    282282                        /*first, retrieve values that are in the f and s-set from the g-set values matrix: */
    283                         localvalues=xNew<double>(this->row_fsize*this->row_ssize);
     283                        localvalues=xNew<IssmDouble>(this->row_fsize*this->row_ssize);
    284284                        for(i=0;i<this->row_fsize;i++){
    285285                                for(j=0;j<this->row_ssize;j++){
    … …  
    291291
    292292                        /*Free ressources:*/
    293                         xDelete<double>(localvalues);
     293                        xDelete<IssmDouble>(localvalues);
    294294                }
    295295        }
    … …  
    304304
    305305        int i,j;
    306         double* localvalues=NULL;
     306        IssmDouble* localvalues=NULL;
    307307
    308308        /*Check that Jff is not NULL*/
    … …  
    317317                if(this->row_fsize){
    318318                        /*first, retrieve values that are in the f-set from the g-set values matrix: */
    319                         localvalues=xNew<double>(this->row_fsize*this->row_fsize);
     319                        localvalues=xNew<IssmDouble>(this->row_fsize*this->row_fsize);
    320320                        for(i=0;i<this->row_fsize;i++){
    321321                                for(j=0;j<this->row_fsize;j++){
    … …  
    327327
    328328                        /*Free ressources:*/
    329                         xDelete<double>(localvalues);
     329                        xDelete<IssmDouble>(localvalues);
    330330                }
    331331
    … …  
    428428        this->dofsymmetrical=Ke->dofsymmetrical;
    429429
    430         this->values=xNew<double>(this->nrows*this->ncols);
    431         memcpy(this->values,Ke->values,this->nrows*this->ncols*sizeof(double));
     430        this->values=xNew<IssmDouble>(this->nrows*this->ncols);
     431        xMemCpy<IssmDouble>(this->values,Ke->values,this->nrows*this->ncols);
    432432
    433433        this->gglobaldoflist=xNew<int>(this->nrows);
    434         memcpy(this->gglobaldoflist,Ke->gglobaldoflist,this->nrows*sizeof(int));
     434        xMemCpy<int>(this->gglobaldoflist,Ke->gglobaldoflist,this->nrows);
    435435
    436436        this->row_fsize=Ke->row_fsize;
    437437        if(this->row_fsize){
    438438                this->row_flocaldoflist=xNew<int>(this->row_fsize);
    439                 memcpy(this->row_flocaldoflist,Ke->row_flocaldoflist,this->row_fsize*sizeof(int));
     439                xMemCpy<int>(this->row_flocaldoflist,Ke->row_flocaldoflist,this->row_fsize);
    440440                this->row_fglobaldoflist=xNew<int>(this->row_fsize);
    441                 memcpy(this->row_fglobaldoflist,Ke->row_fglobaldoflist,this->row_fsize*sizeof(int));
     441                xMemCpy<int>(this->row_fglobaldoflist,Ke->row_fglobaldoflist,this->row_fsize);
    442442        }
    443443        else{
    … …  
    449449        if(this->row_ssize){
    450450                this->row_slocaldoflist=xNew<int>(this->row_ssize);
    451                 memcpy(this->row_slocaldoflist,Ke->row_slocaldoflist,this->row_ssize*sizeof(int));
     451                xMemCpy<int>(this->row_slocaldoflist,Ke->row_slocaldoflist,this->row_ssize);
    452452                this->row_sglobaldoflist=xNew<int>(this->row_ssize);
    453                 memcpy(this->row_sglobaldoflist,Ke->row_sglobaldoflist,this->row_ssize*sizeof(int));
     453                xMemCpy<int>(this->row_sglobaldoflist,Ke->row_sglobaldoflist,this->row_ssize);
    454454        }
    455455        else{
    … …  
    461461        if(this->col_fsize){
    462462                this->col_flocaldoflist=xNew<int>(this->col_fsize);
    463                 memcpy(this->col_flocaldoflist,Ke->col_flocaldoflist,this->col_fsize*sizeof(int));
     463                xMemCpy<int>(this->col_flocaldoflist,Ke->col_flocaldoflist,this->col_fsize);
    464464                this->col_fglobaldoflist=xNew<int>(this->col_fsize);
    465                 memcpy(this->col_fglobaldoflist,Ke->col_fglobaldoflist,this->col_fsize*sizeof(int));
     465                xMemCpy<int>(this->col_fglobaldoflist,Ke->col_fglobaldoflist,this->col_fsize);
    466466        }
    467467        else{
    … …  
    473473        if(this->col_ssize){
    474474                this->col_slocaldoflist=xNew<int>(this->col_ssize);
    475                 memcpy(this->col_slocaldoflist,Ke->col_slocaldoflist,this->col_ssize*sizeof(int));
     475                xMemCpy<int>(this->col_slocaldoflist,Ke->col_slocaldoflist,this->col_ssize);
    476476                this->col_sglobaldoflist=xNew<int>(this->col_ssize);
    477                 memcpy(this->col_sglobaldoflist,Ke->col_sglobaldoflist,this->col_ssize*sizeof(int));
     477                xMemCpy<int>(this->col_sglobaldoflist,Ke->col_sglobaldoflist,this->col_ssize);
    478478        }
    479479        else{
    … …  
    484484/*}}}*/
    485485/*FUNCTION ElementMatrix::SetDiag{{{*/
    486 void ElementMatrix::SetDiag(double scalar){
     486void ElementMatrix::SetDiag(IssmDouble scalar){
    487487
    488488        int i;

  • issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.h

    r12365 r12472  
    2424                int      ncols;
    2525                bool     dofsymmetrical;
    26                 double*  values;
     26                IssmDouble*  values;
    2727
    2828                //gset
    … …  
    6464                void Transpose(void);
    6565                void Init(ElementMatrix* Ke);
    66                 void SetDiag(double scalar);
     66                void SetDiag(IssmDouble scalar);
    6767                /*}}}*/
    6868};

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

    r12446 r12472  
    7575        /*Gset and values*/
    7676        this->gglobaldoflist=xNew<int>(this->nrows);
    77         this->values=xNewZeroInit<double>(this->nrows);
     77        this->values=xNewZeroInit<IssmDouble>(this->nrows);
    7878        for(i=0;i<pe1->nrows;i++){
    7979                this->values[i] += pe1->values[i];
    … …  
    138138
    139139        /*fill values with 0: */
    140         this->values=xNewZeroInit<double>(this->nrows);
     140        this->values=xNewZeroInit<IssmDouble>(this->nrows);
    141141       
    142142        /*g list*/
    … …  
    152152ElementVector::~ElementVector(){
    153153       
    154         xDelete<double>(this->values);
     154        xDelete<IssmDouble>(this->values);
    155155        xDelete<int>(this->gglobaldoflist);
    156156        xDelete<int>(this->flocaldoflist);
    … …  
    164164
    165165        int i;
    166         double* localvalues=NULL;
     166        IssmDouble* localvalues=NULL;
    167167
    168168        /*In debugging mode, check consistency (no NaN, and values not too big)*/
    … …  
    171171        if(this->fsize){
    172172                /*first, retrieve values that are in the f-set from the g-set values vector: */
    173                 localvalues=xNew<double>(this->fsize);
     173                localvalues=xNew<IssmDouble>(this->fsize);
    174174                for(i=0;i<this->fsize;i++){
    175175                        localvalues[i]=this->values[this->flocaldoflist[i]];
    … …  
    179179
    180180                /*Free ressources:*/
    181                 xDelete<double>(localvalues);
     181                xDelete<IssmDouble>(localvalues);
    182182        }
    183183       
    … …  
    188188
    189189        int i;
    190         double* localvalues=NULL;
     190        IssmDouble* localvalues=NULL;
    191191
    192192        if(this->fsize){
    193193                /*first, retrieve values that are in the f-set from the g-set values vector: */
    194                 localvalues=xNew<double>(this->fsize);
     194                localvalues=xNew<IssmDouble>(this->fsize);
    195195                for(i=0;i<this->fsize;i++){
    196196                        localvalues[i]=this->values[this->flocaldoflist[i]];
    … …  
    200200
    201201                /*Free ressources:*/
    202                 xDelete<double>(localvalues);
     202                xDelete<IssmDouble>(localvalues);
    203203        }
    204204
    … …  
    245245        this->nrows =pe->nrows;
    246246
    247         this->values=xNew<double>(this->nrows);
    248         memcpy(this->values,pe->values,this->nrows*sizeof(double));
     247        this->values=xNew<IssmDouble>(this->nrows);
     248        xMemCpy<IssmDouble>(this->values,pe->values,this->nrows);
    249249
    250250        this->gglobaldoflist=xNew<int>(this->nrows);
    251         memcpy(this->gglobaldoflist,pe->gglobaldoflist,this->nrows*sizeof(int));
     251        xMemCpy<int>(this->gglobaldoflist,pe->gglobaldoflist,this->nrows);
    252252
    253253        this->fsize=pe->fsize;
    254254        if(this->fsize){
    255255                this->flocaldoflist=xNew<int>(this->fsize);
    256                 memcpy(this->flocaldoflist,pe->flocaldoflist,this->fsize*sizeof(int));
     256                xMemCpy<int>(this->flocaldoflist,pe->flocaldoflist,this->fsize);
    257257                this->fglobaldoflist=xNew<int>(this->fsize);
    258                 memcpy(this->fglobaldoflist,pe->fglobaldoflist,this->fsize*sizeof(int));
     258                xMemCpy<int>(this->fglobaldoflist,pe->fglobaldoflist,this->fsize);
    259259        }
    260260        else{
    … …  
    265265/*}}}*/
    266266/*FUNCTION ElementVector::SetValue{{{*/
    267 void ElementVector::SetValue(double scalar){
     267void ElementVector::SetValue(IssmDouble scalar){
    268268
    269269        int i;

  • issm/trunk-jpl/src/c/objects/Numerics/ElementVector.h

    r12365 r12472  
    2222       
    2323                int      nrows;
    24                 double*  values;
     24                IssmDouble*  values;
    2525               
    2626                //gset
    … …  
    4545                void CheckConsistency(void);
    4646                void Init(ElementVector* pe);
    47                 void SetValue(double scalar);
     47                void SetValue(IssmDouble scalar);
    4848                /*}}}*/
    4949};

  • issm/trunk-jpl/src/c/objects/Numerics/Matrix.cpp

    r12441 r12472  
    4545        #endif
    4646        #ifdef _HAVE_ADOLC_
    47         this->amatrix=xNew<adouble>(M*N);
    48         #endif
    49 }
    50 /*}}}*/
    51 /*FUNCTION Matrix::Matrix(int M,int N,double sparsity){{{*/
    52 Matrix::Matrix(int M,int N,double sparsity){
     47        this->amatrix=xNew<IssmDouble>(M*N);
     48        #endif
     49}
     50/*}}}*/
     51/*FUNCTION Matrix::Matrix(int M,int N,IssmDouble sparsity){{{*/
     52Matrix::Matrix(int M,int N,IssmDouble sparsity){
    5353
    5454        #ifdef _HAVE_PETSC_
    … …  
    5858        #endif
    5959        #ifdef _HAVE_ADOLC_
    60         this->amatrix=xNew<adouble>;
    61         #endif
    62 }
    63 /*}}}*/
    64 /*FUNCTION Matrix::Matrix(double* serial_mat, int M,int N,double sparsity){{{*/
    65 Matrix::Matrix(double* serial_mat, int M,int N,double sparsity){
     60        this->amatrix=xNew<IssmDouble>;
     61        #endif
     62}
     63/*}}}*/
     64/*FUNCTION Matrix::Matrix(IssmDouble* serial_mat, int M,int N,IssmDouble sparsity){{{*/
     65Matrix::Matrix(IssmDouble* serial_mat, int M,int N,IssmDouble sparsity){
    6666
    6767        #ifdef _HAVE_PETSC_
    … …  
    8585        #endif
    8686        #ifdef _HAVE_ADOLC_
    87         this->amatrix=xNew<adouble>(M*N);
     87        this->amatrix=xNew<IssmDouble>(M*N);
    8888        #endif
    8989}
    … …  
    9898        #endif
    9999        #ifdef _HAVE_ADOLC_
    100         this->amatrix=xNew<adouble>(M*N);
     100        this->amatrix=xNew<IssmDouble>(M*N);
    101101        #endif
    102102}
    … …  
    111111        #endif
    112112        #ifdef _HAVE_ADOLC_
    113         xDelete<adouble>(this->amatrix);
     113        xDelete<IssmDouble>(this->amatrix);
    114114        #endif
    115115}
    … …  
    156156/*}}}*/
    157157/*FUNCTION Matrix::Norm{{{*/
    158 double Matrix::Norm(NormMode norm_type){
    159        
    160         double norm=0;
     158IssmDouble Matrix::Norm(NormMode norm_type){
     159       
     160        IssmDouble norm=0;
    161161        #ifdef _HAVE_PETSC_
    162162                _assert_(this->matrix);
    … …  
    220220/*}}}*/
    221221/*FUNCTION Matrix::ToSerial{{{*/
    222 double* Matrix::ToSerial(void){
    223 
    224         double* output=NULL;
     222IssmDouble* Matrix::ToSerial(void){
     223
     224        IssmDouble* output=NULL;
    225225
    226226        #ifdef _HAVE_PETSC_
    … …  
    233233/*}}}*/
    234234/*FUNCTION Matrix::SetValues{{{*/
    235 void Matrix::SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode){
     235void Matrix::SetValues(int m,int* idxm,int n,int* idxn,IssmDouble* values,InsMode mode){
    236236
    237237        #ifdef _HAVE_PETSC_

  • issm/trunk-jpl/src/c/objects/Numerics/Matrix.h

    r12365 r12472  
    3535                #endif
    3636                #ifdef _HAVE_ADOLC_
    37                 adouble* amatrix;
     37                IssmDouble* amatrix;
    3838                #endif
    3939
    … …  
    4141                Matrix();
    4242                Matrix(int M,int N);
    43                 Matrix(int M,int N,double sparsity);
    44                 Matrix(double* serial_mat,int M,int N,double sparsity);
     43                Matrix(int M,int N,IssmDouble sparsity);
     44                Matrix(IssmDouble* serial_mat,int M,int N,IssmDouble sparsity);
    4545                Matrix(int M,int N,int connectivity,int numberofdofspernode);
    4646                ~Matrix();
    … …  
    4949                void Echo(void);
    5050                void Assemble(void);
    51                 double Norm(NormMode norm_type);
     51                IssmDouble Norm(NormMode norm_type);
    5252                void GetSize(int* pM,int* pN);
    5353                void GetLocalSize(int* pM,int* pN);
    5454                void MatMult(Vector* X,Vector* AX);
    5555                Matrix* Duplicate(void);
    56                 double* ToSerial(void);
    57                 void SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode);
     56                IssmDouble* ToSerial(void);
     57                void SetValues(int m,int* idxm,int n,int* idxn,IssmDouble* values,InsMode mode);
    5858                void Convert(MatrixType type);
    5959                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Numerics/Vector.cpp

    r12441 r12472  
    4444        #endif
    4545        #ifdef _HAVE_ADOLC_
    46         this->avector=xNew<adouble>(pM);
    47         #endif
    48 }
    49 /*}}}*/
    50 /*FUNCTION Vector::Vector(double* serial_vec,int M){{{*/
    51 Vector::Vector(double* serial_vec,int M){
     46        this->avector=xNew<IssmDouble>(pM);
     47        #endif
     48}
     49/*}}}*/
     50/*FUNCTION Vector::Vector(IssmDouble* serial_vec,int M){{{*/
     51Vector::Vector(IssmDouble* serial_vec,int M){
    5252
    5353        #ifdef _HAVE_PETSC_
    … …  
    6565        #endif
    6666        #ifdef _HAVE_ADOLC_
    67                 this->avector=xNew<adouble>(M);
     67                this->avector=xNew<IssmDouble>(M);
    6868        #endif
    6969}
    … …  
    109109        #endif
    110110        #ifdef _HAVE_ADOLC_
    111         xDelete<adouble>(this->avector);
     111        xDelete<IssmDouble>(this->avector);
    112112        #endif
    113113}
    … …  
    146146/*}}}*/
    147147/*FUNCTION Vector::SetValues{{{*/
    148 void Vector::SetValues(int ssize, int* list, double* values, InsMode mode){
     148void Vector::SetValues(int ssize, int* list, IssmDouble* values, InsMode mode){
    149149               
    150150               
    … …  
    159159/*}}}*/
    160160/*FUNCTION Vector::SetValue{{{*/
    161 void Vector::SetValue(int dof, double value, InsMode mode){
     161void Vector::SetValue(int dof, IssmDouble value, InsMode mode){
    162162               
    163163        #ifdef _HAVE_PETSC_
    … …  
    171171/*}}}*/
    172172/*FUNCTION Vector::GetValue{{{*/
    173 void Vector::GetValue(double* pvalue,int dof){
     173void Vector::GetValue(IssmDouble* pvalue,int dof){
    174174               
    175175        #ifdef _HAVE_PETSC_
    … …  
    238238/*}}}*/
    239239/*FUNCTION Vector::Set{{{*/
    240 void Vector::Set(double value){
     240void Vector::Set(IssmDouble value){
    241241       
    242242        #ifdef _HAVE_PETSC_
    … …  
    250250/*}}}*/
    251251/*FUNCTION Vector::AXPY{{{*/
    252 void Vector::AXPY(Vector* X, double a){
     252void Vector::AXPY(Vector* X, IssmDouble a){
    253253       
    254254        #ifdef _HAVE_PETSC_
    … …  
    261261/*}}}*/
    262262/*FUNCTION Vector::AYPX{{{*/
    263 void Vector::AYPX(Vector* X, double a){
     263void Vector::AYPX(Vector* X, IssmDouble a){
    264264       
    265265        #ifdef _HAVE_PETSC_
    … …  
    273273/*}}}*/
    274274/*FUNCTION Vector::ToMPISerial{{{*/
    275 double* Vector::ToMPISerial(void){
    276 
    277         double* vec_serial=NULL;
     275IssmDouble* Vector::ToMPISerial(void){
     276
     277        IssmDouble* vec_serial=NULL;
    278278
    279279        #ifdef _HAVE_PETSC_
    … …  
    299299/*}}}*/
    300300/*FUNCTION Vector::Norm{{{*/
    301 double Vector::Norm(NormMode norm_type){
    302        
    303         double norm=0;
     301IssmDouble Vector::Norm(NormMode norm_type){
     302       
     303        IssmDouble norm=0;
    304304        #ifdef _HAVE_PETSC_
    305305                _assert_(this->vector);
    … …  
    312312/*}}}*/
    313313/*FUNCTION Vector::Scale{{{*/
    314 void Vector::Scale(double scale_factor){
     314void Vector::Scale(IssmDouble scale_factor){
    315315       
    316316        #ifdef _HAVE_PETSC_
    … …  
    323323/*}}}*/
    324324/*FUNCTION Vector::Dot{{{*/
    325 double Vector::Dot(Vector* vector){
    326 
    327         double dot;
     325IssmDouble Vector::Dot(Vector* vector){
     326
     327        IssmDouble dot;
    328328        #ifdef _HAVE_PETSC_
    329329                _assert_(this->vector);

  • issm/trunk-jpl/src/c/objects/Numerics/Vector.h

    r12365 r12472  
    3333                #endif
    3434                #ifdef _HAVE_ADOLC_
    35                         adouble* avector;
     35                        IssmDouble* avector;
    3636                #endif
    3737
    … …  
    3939                Vector();
    4040                Vector(int M,bool fromlocalsize=false);
    41                 Vector(double* serial_vec,int pM);
     41                Vector(IssmDouble* serial_vec,int pM);
    4242                #ifdef _HAVE_PETSC_
    4343                Vector(Vec petsc_vec);
    … …  
    5050                /*Vector specific routines {{{*/
    5151                void Echo(void);
    52                 void    AXPY(Vector *X, double a);
    53                 void    AYPX(Vector *X, double a);
     52                void    AXPY(Vector *X, IssmDouble a);
     53                void    AYPX(Vector *X, IssmDouble a);
    5454                void    Assemble(void);
    5555                void    Copy(Vector *to);
    56                 double  Dot(Vector *vector);
     56                IssmDouble  Dot(Vector *vector);
    5757                Vector *Duplicate(void);
    58                 void    GetValue(double *pvalue, int dof);
     58                void    GetValue(IssmDouble *pvalue, int dof);
    5959                void    GetSize(int *pM);
    6060                void    GetLocalSize(int *pM);
    6161                bool    IsEmpty(void);
    62                 double  Norm(NormMode norm_type);
     62                IssmDouble  Norm(NormMode norm_type);
    6363                void    PointwiseDivide(Vector  *x,Vector*y);
    64                 void    Scale(double scale_factor);
    65                 void    Set(double value);
    66                 void    SetValues(int ssize, int *list, double*values, InsMode mode);
    67                 void    SetValue(int dof, double value, InsMode mode);
    68                 double *ToMPISerial(void);
     64                void    Scale(IssmDouble scale_factor);
     65                void    Set(IssmDouble value);
     66                void    SetValues(int ssize, int *list, IssmDouble*values, InsMode mode);
     67                void    SetValue(int dof, IssmDouble value, InsMode mode);
     68                IssmDouble *ToMPISerial(void);
    6969                /*}}}*/
    7070};
Note: See TracChangeset for help on using the changeset viewer.