Index: /issm/trunk-jpl/src/c/objects/Loads/Friction.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Friction.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Friction.cpp (revision 12472) @@ -34,5 +34,5 @@ this->inputs=inputs_in; this->element_type=xNew(strlen(element_type_in)+1); - memcpy(this->element_type,element_type_in,(strlen(element_type_in)+1)*sizeof(char)); + xMemCpy(this->element_type,element_type_in,(strlen(element_type_in)+1)); this->matpar=matpar_in; Index: /issm/trunk-jpl/src/c/objects/Loads/Icefront.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Icefront.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Icefront.cpp (revision 12472) @@ -279,5 +279,5 @@ /*}}}*/ /*FUNCTION Icefront::PenaltyCreateKMatrix {{{*/ -void Icefront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, double kmax){ +void Icefront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, IssmDouble kmax){ /*do nothing: */ return; @@ -285,5 +285,5 @@ /*}}}*/ /*FUNCTION Icefront::PenaltyCreatePVector{{{*/ -void Icefront::PenaltyCreatePVector(Vector* pf,double kmax){ +void Icefront::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){ /*do nothing: */ return; @@ -291,5 +291,5 @@ /*}}}*/ /*FUNCTION Icefront::PenaltyCreateJacobianMatrix{{{*/ -void Icefront::PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){ +void Icefront::PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){ this->PenaltyCreateKMatrix(Jff,NULL,kmax); } @@ -303,6 +303,6 @@ /*Update virtual functions definitions:*/ -/*FUNCTION Icefront::InputUpdateFromVector(double* vector, int name, int type) {{{*/ -void Icefront::InputUpdateFromVector(double* vector, int name, int type){ +/*FUNCTION Icefront::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/ +void Icefront::InputUpdateFromVector(IssmDouble* vector, int name, int type){ /*Nothing updated yet*/ } @@ -318,11 +318,11 @@ } /*}}}*/ -/*FUNCTION Icefront::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols, int name, int type) {{{*/ -void Icefront::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols, int name, int type){ - /*Nothing updated yet*/ -} -/*}}}*/ -/*FUNCTION Icefront::InputUpdateFromVectorDakota(double* vector, int name, int type) {{{*/ -void Icefront::InputUpdateFromVectorDakota(double* vector, int name, int type){ +/*FUNCTION Icefront::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type) {{{*/ +void Icefront::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type){ + /*Nothing updated yet*/ +} +/*}}}*/ +/*FUNCTION Icefront::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type) {{{*/ +void Icefront::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){ /*Nothing updated yet*/ } @@ -338,6 +338,6 @@ } /*}}}*/ -/*FUNCTION Icefront::InputUpdateFromConstant(double constant, int name) {{{*/ -void Icefront::InputUpdateFromConstant(double constant, int name){ +/*FUNCTION Icefront::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/ +void Icefront::InputUpdateFromConstant(IssmDouble constant, int name){ /*Nothing updated yet*/ } @@ -354,5 +354,5 @@ /*}}}*/ /*FUNCTION Icefront::InputUpdateFromSolution{{{*/ -void Icefront::InputUpdateFromSolution(double* solution){ +void Icefront::InputUpdateFromSolution(IssmDouble* solution){ /*Nothing updated yet*/ } @@ -392,10 +392,10 @@ /*Intermediary*/ int ig,index1,index2,fill; - double Jdet; - double thickness,bed,pressure,ice_pressure,rho_water,rho_ice,gravity; - double water_pressure,air_pressure,surface_under_water,base_under_water; - double xyz_list[numnodes][3]; - double normal[2]; - double L[2]; + IssmDouble Jdet; + IssmDouble thickness,bed,pressure,ice_pressure,rho_water,rho_ice,gravity; + IssmDouble water_pressure,air_pressure,surface_under_water,base_under_water; + IssmDouble xyz_list[numnodes][3]; + IssmDouble normal[2]; + IssmDouble L[2]; GaussTria *gauss; @@ -513,10 +513,10 @@ int i,j,ig,index1,index2,index3,index4; int fill; - double surface,pressure,ice_pressure,rho_water,rho_ice,gravity; - double water_pressure,air_pressure; - double Jdet,z_g; - double xyz_list[NUMVERTICESQUA][3]; - double normal[3]; - double l1l4[4]; + IssmDouble surface,pressure,ice_pressure,rho_water,rho_ice,gravity; + IssmDouble water_pressure,air_pressure; + IssmDouble Jdet,z_g; + IssmDouble xyz_list[NUMVERTICESQUA][3]; + IssmDouble normal[3]; + IssmDouble l1l4[4]; GaussPenta *gauss = NULL; @@ -543,6 +543,6 @@ /* Start looping on the number of gaussian points: */ - double zmax=xyz_list[0][2]; for(i=1;izmax) zmax=xyz_list[i][2]; - double zmin=xyz_list[0][2]; for(i=1;izmax) zmax=xyz_list[i][2]; + IssmDouble zmin=xyz_list[0][2]; for(i=1;i0 && zmin<0) gauss=new GaussPenta(index1,index2,index3,index4,3,10); //refined in vertical because of the sea level discontinuity else gauss=new GaussPenta(index1,index2,index3,index4,3,3); @@ -590,10 +590,10 @@ int i,j,ig,index1,index2,index3,index4; int fill; - double pressure,rho_water,gravity; - double water_pressure,air_pressure; - double Jdet,z_g; - double xyz_list[NUMVERTICESQUA][3]; - double normal[3]; - double l1l4[4]; + IssmDouble pressure,rho_water,gravity; + IssmDouble water_pressure,air_pressure; + IssmDouble Jdet,z_g; + IssmDouble xyz_list[NUMVERTICESQUA][3]; + IssmDouble normal[3]; + IssmDouble l1l4[4]; GaussPenta *gauss = NULL; @@ -618,6 +618,6 @@ /* Start looping on the number of gaussian points: */ - double zmax=xyz_list[0][2]; for(i=1;izmax) zmax=xyz_list[i][2]; - double zmin=xyz_list[0][2]; for(i=1;izmax) zmax=xyz_list[i][2]; + IssmDouble zmin=xyz_list[0][2]; for(i=1;i0 && zmin<0) gauss=new GaussPenta(index1,index2,index3,index4,3,30); //refined in vertical because of the sea level discontinuity else gauss=new GaussPenta(index1,index2,index3,index4,3,3); @@ -705,10 +705,10 @@ /*}}}*/ /*FUNCTION Icefront::GetSegmentNormal {{{*/ -void Icefront:: GetSegmentNormal(double* normal,double xyz_list[4][3]){ +void Icefront:: GetSegmentNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]){ /*Build unit outward pointing vector*/ const int numnodes=NUMVERTICESSEG; - double vector[2]; - double norm; + IssmDouble vector[2]; + IssmDouble norm; vector[0]=xyz_list[1][0] - xyz_list[0][0]; @@ -722,10 +722,10 @@ /*}}}*/ /*FUNCTION Icefront::GetQuadNormal {{{*/ -void Icefront:: GetQuadNormal(double* normal,double xyz_list[4][3]){ +void Icefront:: GetQuadNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]){ /*Build unit outward pointing vector*/ - double AB[3]; - double AC[3]; - double norm; + IssmDouble AB[3]; + IssmDouble AC[3]; + IssmDouble norm; AB[0]=xyz_list[1][0] - xyz_list[0][0]; Index: /issm/trunk-jpl/src/c/objects/Loads/Icefront.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Icefront.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Icefront.h (revision 12472) @@ -53,15 +53,15 @@ /*}}}*/ /*Update virtual functions definitions: {{{*/ - void InputUpdateFromVector(double* vector, int name, int type); + void InputUpdateFromVector(IssmDouble* vector, int name, int type); void InputUpdateFromVector(int* vector, int name, int type); void InputUpdateFromVector(bool* vector, int name, int type); - void InputUpdateFromMatrixDakota(double* matrix,int ncols,int nrows, int name, int type); - void InputUpdateFromVectorDakota(double* vector, int name, int type); + void InputUpdateFromMatrixDakota(IssmDouble* matrix,int ncols,int nrows, int name, int type); + void InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type); void InputUpdateFromVectorDakota(int* vector, int name, int type); void InputUpdateFromVectorDakota(bool* vector, int name, int type); - void InputUpdateFromConstant(double constant, int name); + void InputUpdateFromConstant(IssmDouble constant, int name); void InputUpdateFromConstant(int constant, int name); void InputUpdateFromConstant(bool constant, int name); - void InputUpdateFromSolution(double* solution); + void InputUpdateFromSolution(IssmDouble* solution); void InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");}; /*}}}*/ @@ -72,13 +72,13 @@ void CreatePVector(Vector* pf); void CreateJacobianMatrix(Matrix* Jff); - void PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax); - void PenaltyCreatePVector(Vector* pf, double kmax); - void PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax); + void PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, IssmDouble kmax); + void PenaltyCreatePVector(Vector* pf, IssmDouble kmax); + void PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax); bool InAnalysis(int analysis_type); /*}}}*/ /*Load management: {{{*/ void GetDofList(int** pdoflist,int approximation_enum,int setenum); - void GetSegmentNormal(double* normal,double xyz_list[2][3]); - void GetQuadNormal(double* normal,double xyz_list[4][3]); + void GetSegmentNormal(IssmDouble* normal,IssmDouble xyz_list[2][3]); + void GetQuadNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]); #ifdef _HAVE_CONTROL_ ElementVector* CreatePVectorAdjointHoriz(void); Index: /issm/trunk-jpl/src/c/objects/Loads/Load.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Load.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Load.h (revision 12472) @@ -32,7 +32,7 @@ virtual void CreatePVector(Vector* pf)=0; virtual void CreateJacobianMatrix(Matrix* Jff)=0; - virtual void PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax)=0; - virtual void PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, double kmax)=0; - virtual void PenaltyCreatePVector(Vector* pf, double kmax)=0; + virtual void PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax)=0; + virtual void PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, IssmDouble kmax)=0; + virtual void PenaltyCreatePVector(Vector* pf, IssmDouble kmax)=0; virtual bool InAnalysis(int analysis_type)=0; /*}}}*/ Index: /issm/trunk-jpl/src/c/objects/Loads/Numericalflux.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Numericalflux.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Numericalflux.cpp (revision 12472) @@ -60,5 +60,5 @@ /*First, see wether this is an internal or boundary edge (if e2=NaN)*/ - if (isnan((double)iomodel->Data(MeshEdgesEnum)[4*i+3])){ //edges are [node1 node2 elem1 elem2] + if (isnan((IssmDouble)iomodel->Data(MeshEdgesEnum)[4*i+3])){ //edges are [node1 node2 elem1 elem2] /* Boundary edge, only one element */ e1=(int)iomodel->Data(MeshEdgesEnum)[4*i+2]; @@ -312,5 +312,5 @@ /*}}}*/ /*FUNCTION Numericalflux::PenaltyCreateKMatrix {{{*/ -void Numericalflux::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){ +void Numericalflux::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,IssmDouble kmax){ /*No stiffness loads applied, do nothing: */ @@ -320,5 +320,5 @@ /*}}}*/ /*FUNCTION Numericalflux::PenaltyCreatePVector{{{*/ -void Numericalflux::PenaltyCreatePVector(Vector* pf,double kmax){ +void Numericalflux::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){ /*No penalty loads applied, do nothing: */ @@ -359,11 +359,11 @@ /* Intermediaries*/ int i,j,ig,index1,index2; - double DL1,DL2,Jdet,dt,vx,vy,UdotN; - double xyz_list[NUMVERTICES_INTERNAL][3]; - double normal[2]; - double B[numdof]; - double Bprime[numdof]; - double Ke_g1[numdof][numdof]; - double Ke_g2[numdof][numdof]; + IssmDouble DL1,DL2,Jdet,dt,vx,vy,UdotN; + IssmDouble xyz_list[NUMVERTICES_INTERNAL][3]; + IssmDouble normal[2]; + IssmDouble B[numdof]; + IssmDouble Bprime[numdof]; + IssmDouble Ke_g1[numdof][numdof]; + IssmDouble Ke_g2[numdof][numdof]; GaussTria *gauss; @@ -424,9 +424,9 @@ /* Intermediaries*/ int i,j,ig,index1,index2; - double DL,Jdet,dt,vx,vy,mean_vx,mean_vy,UdotN; - double xyz_list[NUMVERTICES_BOUNDARY][3]; - double normal[2]; - double L[numdof]; - double Ke_g[numdof][numdof]; + IssmDouble DL,Jdet,dt,vx,vy,mean_vx,mean_vy,UdotN; + IssmDouble xyz_list[NUMVERTICES_BOUNDARY][3]; + IssmDouble normal[2]; + IssmDouble L[numdof]; + IssmDouble Ke_g[numdof][numdof]; GaussTria *gauss; @@ -512,11 +512,11 @@ /* Intermediaries*/ int i,j,ig,index1,index2; - double DL1,DL2,Jdet,vx,vy,UdotN; - double xyz_list[NUMVERTICES_INTERNAL][3]; - double normal[2]; - double B[numdof]; - double Bprime[numdof]; - double Ke_g1[numdof][numdof]; - double Ke_g2[numdof][numdof]; + IssmDouble DL1,DL2,Jdet,vx,vy,UdotN; + IssmDouble xyz_list[NUMVERTICES_INTERNAL][3]; + IssmDouble normal[2]; + IssmDouble B[numdof]; + IssmDouble Bprime[numdof]; + IssmDouble Ke_g1[numdof][numdof]; + IssmDouble Ke_g2[numdof][numdof]; GaussTria *gauss; @@ -576,9 +576,9 @@ /* Intermediaries*/ int i,j,ig,index1,index2; - double DL,Jdet,vx,vy,mean_vx,mean_vy,UdotN; - double xyz_list[NUMVERTICES_BOUNDARY][3]; - double normal[2]; - double L[numdof]; - double Ke_g[numdof][numdof]; + IssmDouble DL,Jdet,vx,vy,mean_vx,mean_vy,UdotN; + IssmDouble xyz_list[NUMVERTICES_BOUNDARY][3]; + IssmDouble normal[2]; + IssmDouble L[numdof]; + IssmDouble Ke_g[numdof][numdof]; GaussTria *gauss; @@ -703,8 +703,8 @@ /* Intermediaries*/ int i,j,ig,index1,index2; - double DL,Jdet,dt,vx,vy,mean_vx,mean_vy,UdotN,thickness; - double xyz_list[NUMVERTICES_BOUNDARY][3]; - double normal[2]; - double L[numdof]; + IssmDouble DL,Jdet,dt,vx,vy,mean_vx,mean_vy,UdotN,thickness; + IssmDouble xyz_list[NUMVERTICES_BOUNDARY][3]; + IssmDouble normal[2]; + IssmDouble L[numdof]; GaussTria *gauss; @@ -797,8 +797,8 @@ /* Intermediaries*/ int i,j,ig,index1,index2; - double DL,Jdet,vx,vy,mean_vx,mean_vy,UdotN,thickness; - double xyz_list[NUMVERTICES_BOUNDARY][3]; - double normal[2]; - double L[numdof]; + IssmDouble DL,Jdet,vx,vy,mean_vx,mean_vy,UdotN,thickness; + IssmDouble xyz_list[NUMVERTICES_BOUNDARY][3]; + IssmDouble normal[2]; + IssmDouble L[numdof]; GaussTria *gauss; @@ -864,9 +864,9 @@ /*}}}*/ /*FUNCTION Numericalflux::GetNormal {{{*/ -void Numericalflux:: GetNormal(double* normal,double xyz_list[4][3]){ +void Numericalflux:: GetNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]){ /*Build unit outward pointing vector*/ - double vector[2]; - double norm; + IssmDouble vector[2]; + IssmDouble norm; vector[0]=xyz_list[1][0] - xyz_list[0][0]; Index: /issm/trunk-jpl/src/c/objects/Loads/Numericalflux.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Numericalflux.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Numericalflux.h (revision 12472) @@ -49,15 +49,15 @@ /*}}}*/ /*Update virtual functions resolution: {{{*/ - void InputUpdateFromVector(double* vector, int name, int type){/*Do nothing*/} + void InputUpdateFromVector(IssmDouble* vector, int name, int type){/*Do nothing*/} void InputUpdateFromVector(int* vector, int name, int type){_error_("Not implemented yet!");} void InputUpdateFromVector(bool* vector, int name, int type){_error_("Not implemented yet!");} - void InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols, int name, int type){/*Do nothing*/} - void InputUpdateFromVectorDakota(double* vector, int name, int type){/*Do nothing*/} + void InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type){/*Do nothing*/} + void InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){/*Do nothing*/} void InputUpdateFromVectorDakota(int* vector, int name, int type){_error_("Not implemented yet!");} void InputUpdateFromVectorDakota(bool* vector, int name, int type){_error_("Not implemented yet!");} - void InputUpdateFromConstant(double constant, int name){/*Do nothing*/}; + void InputUpdateFromConstant(IssmDouble constant, int name){/*Do nothing*/}; void InputUpdateFromConstant(int constant, int name){/*Do nothing*/}; void InputUpdateFromConstant(bool constant, int name){_error_("Not implemented yet!");} - void InputUpdateFromSolution(double* solution){_error_("Not implemented yet!");} + void InputUpdateFromSolution(IssmDouble* solution){_error_("Not implemented yet!");} void InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");}; /*}}}*/ @@ -68,11 +68,11 @@ void CreatePVector(Vector* pf); void CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");}; - void PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");}; - void PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax); - void PenaltyCreatePVector(Vector* pf, double kmax); + void PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){_error_("Not implemented yet");}; + void PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, IssmDouble kmax); + void PenaltyCreatePVector(Vector* pf, IssmDouble kmax); bool InAnalysis(int analysis_type); /*}}}*/ /*Numericalflux management:{{{*/ - void GetNormal(double* normal,double xyz_list[4][3]); + void GetNormal(IssmDouble* normal,IssmDouble xyz_list[4][3]); ElementMatrix* CreateKMatrixPrognostic(void); ElementMatrix* CreateKMatrixPrognosticInternal(void); Index: /issm/trunk-jpl/src/c/objects/Loads/Pengrid.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Pengrid.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Pengrid.cpp (revision 12472) @@ -217,5 +217,5 @@ /*}}}*/ /*FUNCTION Pengrid::PenaltyCreateMatrix {{{*/ -void Pengrid::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){ +void Pengrid::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,IssmDouble kmax){ /*Retrieve parameters: */ @@ -250,5 +250,5 @@ /*}}}*/ /*FUNCTION Pengrid::PenaltyCreatePVector {{{*/ -void Pengrid::PenaltyCreatePVector(Vector* pf,double kmax){ +void Pengrid::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){ /*Retrieve parameters: */ @@ -289,6 +289,6 @@ /*Update virtual functions definitions:*/ -/*FUNCTION Pengrid::InputUpdateFromVector(double* vector, int name, int type) {{{*/ -void Pengrid::InputUpdateFromVector(double* vector, int name, int type){ +/*FUNCTION Pengrid::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/ +void Pengrid::InputUpdateFromVector(IssmDouble* vector, int name, int type){ /*Nothing updated yet*/ } @@ -304,11 +304,11 @@ } /*}}}*/ -/*FUNCTION Pengrid::InputUpdateFromMatrixDakota(double* vector, int nrows, int ncols, int name, int type) {{{*/ -void Pengrid::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols, int name, int type){ - /*Nothing updated yet*/ -} -/*}}}*/ -/*FUNCTION Pengrid::InputUpdateFromVectorDakota(double* vector, int name, int type) {{{*/ -void Pengrid::InputUpdateFromVectorDakota(double* vector, int name, int type){ +/*FUNCTION Pengrid::InputUpdateFromMatrixDakota(IssmDouble* vector, int nrows, int ncols, int name, int type) {{{*/ +void Pengrid::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type){ + /*Nothing updated yet*/ +} +/*}}}*/ +/*FUNCTION Pengrid::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type) {{{*/ +void Pengrid::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){ /*Nothing updated yet*/ } @@ -324,6 +324,6 @@ } /*}}}*/ -/*FUNCTION Pengrid::InputUpdateFromConstant(double constant, int name) {{{*/ -void Pengrid::InputUpdateFromConstant(double constant, int name){ +/*FUNCTION Pengrid::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/ +void Pengrid::InputUpdateFromConstant(IssmDouble constant, int name){ switch(name){ @@ -353,5 +353,5 @@ /*}}}*/ /*FUNCTION Pengrid::InputUpdateFromSolution{{{*/ -void Pengrid::InputUpdateFromSolution(double* solution){ +void Pengrid::InputUpdateFromSolution(IssmDouble* solution){ /*Nothing updated yet*/ } @@ -391,7 +391,7 @@ int found=0; const int numnodes=1; - double pressure; - double temperature; - double t_pmp; + IssmDouble pressure; + IssmDouble temperature; + IssmDouble t_pmp; int new_active; int unstable=0; @@ -455,9 +455,9 @@ #ifdef _HAVE_DIAGNOSTIC_ /*FUNCTION Pengrid::PenaltyCreateKMatrixDiagnosticStokes {{{*/ -ElementMatrix* Pengrid::PenaltyCreateKMatrixDiagnosticStokes(double kmax){ +ElementMatrix* Pengrid::PenaltyCreateKMatrixDiagnosticStokes(IssmDouble kmax){ const int numdof = NUMVERTICES *NDOF4; - double slope[2]; - double penalty_offset; + IssmDouble slope[2]; + IssmDouble penalty_offset; int approximation; @@ -475,7 +475,7 @@ /*Create elementary matrix: add penalty to constrain wb (wb=ub*db/dx+vb*db/dy)*/ - Ke->values[2*NDOF4+0]=-slope[0]*kmax*pow((double)10.0,penalty_offset); - Ke->values[2*NDOF4+1]=-slope[1]*kmax*pow((double)10.0,penalty_offset); - Ke->values[2*NDOF4+2]= kmax*pow((double)10,penalty_offset); + Ke->values[2*NDOF4+0]=-slope[0]*kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[2*NDOF4+1]=-slope[1]*kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[2*NDOF4+2]= kmax*pow((IssmDouble)10,penalty_offset); /*Transform Coordinate System*/ @@ -489,9 +489,9 @@ #ifdef _HAVE_THERMAL_ /*FUNCTION Pengrid::PenaltyCreateKMatrixMelting {{{*/ -ElementMatrix* Pengrid::PenaltyCreateKMatrixMelting(double kmax){ +ElementMatrix* Pengrid::PenaltyCreateKMatrixMelting(IssmDouble kmax){ const int numdof=NUMVERTICES*NDOF1; - double pressure,temperature,t_pmp; - double penalty_factor; + IssmDouble pressure,temperature,t_pmp; + IssmDouble penalty_factor; Penta* penta=(Penta*)element; @@ -511,5 +511,5 @@ /*Add penalty load*/ if (temperaturevalues[0]=kmax*pow((double)10,penalty_factor); + Ke->values[0]=kmax*pow((IssmDouble)10,penalty_factor); } @@ -519,8 +519,8 @@ /*}}}*/ /*FUNCTION Pengrid::PenaltyCreateKMatrixThermal {{{*/ -ElementMatrix* Pengrid::PenaltyCreateKMatrixThermal(double kmax){ +ElementMatrix* Pengrid::PenaltyCreateKMatrixThermal(IssmDouble kmax){ const int numdof=NUMVERTICES*NDOF1; - double penalty_factor; + IssmDouble penalty_factor; /*Initialize Element matrix and return if necessary*/ @@ -531,5 +531,5 @@ parameters->FindParam(&penalty_factor,ThermalPenaltyFactorEnum); - Ke->values[0]=kmax*pow((double)10,penalty_factor); + Ke->values[0]=kmax*pow((IssmDouble)10,penalty_factor); /*Clean up and return*/ @@ -538,12 +538,12 @@ /*}}}*/ /*FUNCTION Pengrid::PenaltyCreatePVectorMelting {{{*/ -ElementVector* Pengrid::PenaltyCreatePVectorMelting(double kmax){ +ElementVector* Pengrid::PenaltyCreatePVectorMelting(IssmDouble kmax){ const int numdof=NUMVERTICES*NDOF1; - double pressure; - double temperature; - double melting_offset; - double t_pmp; - double dt,penalty_factor; + IssmDouble pressure; + IssmDouble temperature; + IssmDouble melting_offset; + IssmDouble t_pmp; + IssmDouble dt,penalty_factor; /*recover pointers: */ @@ -572,6 +572,6 @@ } else{ - if (dt) pe->values[0]=melting_offset*pow((double)10,penalty_factor)*(temperature-t_pmp)/dt; - else pe->values[0]=melting_offset*pow((double)10,penalty_factor)*(temperature-t_pmp); + if (dt) pe->values[0]=melting_offset*pow((IssmDouble)10,penalty_factor)*(temperature-t_pmp)/dt; + else pe->values[0]=melting_offset*pow((IssmDouble)10,penalty_factor)*(temperature-t_pmp); } @@ -581,10 +581,10 @@ /*}}}*/ /*FUNCTION Pengrid::PenaltyCreatePVectorThermal {{{*/ -ElementVector* Pengrid::PenaltyCreatePVectorThermal(double kmax){ +ElementVector* Pengrid::PenaltyCreatePVectorThermal(IssmDouble kmax){ const int numdof=NUMVERTICES*NDOF1; - double pressure; - double t_pmp; - double penalty_factor; + IssmDouble pressure; + IssmDouble t_pmp; + IssmDouble penalty_factor; Penta* penta=(Penta*)element; @@ -601,5 +601,5 @@ t_pmp=matpar->GetMeltingPoint()-matpar->GetBeta()*pressure; - pe->values[0]=kmax*pow((double)10,penalty_factor)*t_pmp; + pe->values[0]=kmax*pow((IssmDouble)10,penalty_factor)*t_pmp; /*Clean up and return*/ @@ -615,5 +615,5 @@ /*}}}*/ /*FUNCTION Pengrid::UpdateInputs {{{*/ -void Pengrid::UpdateInputs(double* solution){ +void Pengrid::UpdateInputs(IssmDouble* solution){ _error_("not supported yet!"); } Index: /issm/trunk-jpl/src/c/objects/Loads/Pengrid.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Pengrid.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Pengrid.h (revision 12472) @@ -54,15 +54,15 @@ /*}}}*/ /*Update virtual functions resolution: {{{*/ - void InputUpdateFromVector(double* vector, int name, int type); + void InputUpdateFromVector(IssmDouble* vector, int name, int type); void InputUpdateFromVector(int* vector, int name, int type); void InputUpdateFromVector(bool* vector, int name, int type); - void InputUpdateFromMatrixDakota(double* matrix ,int nrows, int ncols, int name, int type); - void InputUpdateFromVectorDakota(double* vector, int name, int type); + void InputUpdateFromMatrixDakota(IssmDouble* matrix ,int nrows, int ncols, int name, int type); + void InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type); void InputUpdateFromVectorDakota(int* vector, int name, int type); void InputUpdateFromVectorDakota(bool* vector, int name, int type); - void InputUpdateFromConstant(double constant, int name); + void InputUpdateFromConstant(IssmDouble constant, int name); void InputUpdateFromConstant(int constant, int name); void InputUpdateFromConstant(bool constant, int name); - void InputUpdateFromSolution(double* solution); + void InputUpdateFromSolution(IssmDouble* solution); void InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");}; /*}}}*/ @@ -73,22 +73,22 @@ void CreatePVector(Vector* pf); void CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");}; - void PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");}; - void PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax); - void PenaltyCreatePVector(Vector* pf, double kmax); + void PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){_error_("Not implemented yet");}; + void PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, IssmDouble kmax); + void PenaltyCreatePVector(Vector* pf, IssmDouble kmax); bool InAnalysis(int analysis_type); /*}}}*/ /*Pengrid management {{{*/ #ifdef _HAVE_DIAGNOSTIC_ - ElementMatrix* PenaltyCreateKMatrixDiagnosticStokes(double kmax); + ElementMatrix* PenaltyCreateKMatrixDiagnosticStokes(IssmDouble kmax); #endif #ifdef _HAVE_THERMAL_ - ElementMatrix* PenaltyCreateKMatrixThermal(double kmax); - ElementMatrix* PenaltyCreateKMatrixMelting(double kmax); - ElementVector* PenaltyCreatePVectorThermal(double kmax); - ElementVector* PenaltyCreatePVectorMelting(double kmax); + ElementMatrix* PenaltyCreateKMatrixThermal(IssmDouble kmax); + ElementMatrix* PenaltyCreateKMatrixMelting(IssmDouble kmax); + ElementVector* PenaltyCreatePVectorThermal(IssmDouble kmax); + ElementVector* PenaltyCreatePVectorMelting(IssmDouble kmax); #endif void ConstraintActivate(int* punstable); void ConstraintActivateThermal(int* punstable); - void UpdateInputs(double* solution); + void UpdateInputs(IssmDouble* solution); void ResetConstraint(void); /*}}}*/ Index: /issm/trunk-jpl/src/c/objects/Loads/Penpair.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Penpair.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Penpair.cpp (revision 12472) @@ -157,5 +157,5 @@ /*}}}*/ /*FUNCTION Penpair::PenaltyCreateKMatrix {{{*/ -void Penpair::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){ +void Penpair::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,IssmDouble kmax){ /*Retrieve parameters: */ @@ -183,5 +183,5 @@ /*}}}*/ /*FUNCTION Penpair::PenaltyCreatePVector {{{*/ -void Penpair::PenaltyCreatePVector(Vector* pf,double kmax){ +void Penpair::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){ /*No loads applied, do nothing: */ return; @@ -189,5 +189,5 @@ /*}}}*/ /*FUNCTION Penpair::PenaltyCreateJacobianMatrix{{{*/ -void Penpair::PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){ +void Penpair::PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){ this->PenaltyCreateKMatrix(Jff,NULL,kmax); } @@ -201,6 +201,6 @@ /*Update virtual functions definitions:*/ -/*FUNCTION Penpair::InputUpdateFromConstant(double constant, int name) {{{*/ -void Penpair::InputUpdateFromConstant(double constant, int name){ +/*FUNCTION Penpair::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/ +void Penpair::InputUpdateFromConstant(IssmDouble constant, int name){ /*Nothing updated yet*/ } @@ -216,6 +216,6 @@ } /*}}}*/ -/*FUNCTION Penpair::InputUpdateFromVector(double* vector, int name, int type) {{{*/ -void Penpair::InputUpdateFromVector(double* vector, int name, int type){ +/*FUNCTION Penpair::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/ +void Penpair::InputUpdateFromVector(IssmDouble* vector, int name, int type){ /*Nothing updated yet*/ } @@ -234,5 +234,5 @@ /*Penpair management:*/ /*FUNCTION Penpair::PenaltyCreateKMatrixDiagnosticHoriz{{{*/ -ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticHoriz(double kmax){ +ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticHoriz(IssmDouble kmax){ int approximation0=nodes[0]->GetApproximation(); @@ -269,8 +269,8 @@ /*}}}*/ /*FUNCTION Penpair::PenaltyCreateKMatrixDiagnosticMacAyealPattyn {{{*/ -ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticMacAyealPattyn(double kmax){ +ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticMacAyealPattyn(IssmDouble kmax){ const int numdof=NUMVERTICES*NDOF2; - double penalty_offset; + IssmDouble penalty_offset; /*Initialize Element vector and return if necessary*/ @@ -281,13 +281,13 @@ //Create elementary matrix: add penalty to - Ke->values[0*numdof+0]=+kmax*pow((double)10.0,penalty_offset); - Ke->values[0*numdof+2]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[2*numdof+0]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[2*numdof+2]=+kmax*pow((double)10.0,penalty_offset); - - Ke->values[1*numdof+1]=+kmax*pow((double)10.0,penalty_offset); - Ke->values[1*numdof+3]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[3*numdof+1]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[3*numdof+3]=+kmax*pow((double)10.0,penalty_offset); + Ke->values[0*numdof+0]=+kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[0*numdof+2]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[2*numdof+0]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[2*numdof+2]=+kmax*pow((IssmDouble)10.0,penalty_offset); + + Ke->values[1*numdof+1]=+kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[1*numdof+3]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[3*numdof+1]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[3*numdof+3]=+kmax*pow((IssmDouble)10.0,penalty_offset); /*Clean up and return*/ @@ -296,8 +296,8 @@ /*}}}*/ /*FUNCTION Penpair::PenaltyCreateKMatrixDiagnosticStokes {{{*/ -ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticStokes(double kmax){ +ElementMatrix* Penpair::PenaltyCreateKMatrixDiagnosticStokes(IssmDouble kmax){ const int numdof=NUMVERTICES*NDOF4; - double penalty_offset; + IssmDouble penalty_offset; /*Initialize Element vector and return if necessary*/ @@ -308,23 +308,23 @@ //Create elementary matrix: add penalty to - Ke->values[0*numdof+0]=+kmax*pow((double)10.0,penalty_offset); - Ke->values[0*numdof+4]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[4*numdof+0]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[4*numdof+4]=+kmax*pow((double)10.0,penalty_offset); - - Ke->values[1*numdof+1]=+kmax*pow((double)10.0,penalty_offset); - Ke->values[1*numdof+5]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[5*numdof+1]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[5*numdof+5]=+kmax*pow((double)10.0,penalty_offset); - - Ke->values[2*numdof+2]=+kmax*pow((double)10.0,penalty_offset); - Ke->values[2*numdof+6]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[6*numdof+2]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[6*numdof+6]=+kmax*pow((double)10.0,penalty_offset); - - Ke->values[3*numdof+3]=+kmax*pow((double)10.0,penalty_offset); - Ke->values[3*numdof+7]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[7*numdof+3]=-kmax*pow((double)10.0,penalty_offset); - Ke->values[7*numdof+7]=+kmax*pow((double)10.0,penalty_offset); + Ke->values[0*numdof+0]=+kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[0*numdof+4]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[4*numdof+0]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[4*numdof+4]=+kmax*pow((IssmDouble)10.0,penalty_offset); + + Ke->values[1*numdof+1]=+kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[1*numdof+5]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[5*numdof+1]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[5*numdof+5]=+kmax*pow((IssmDouble)10.0,penalty_offset); + + Ke->values[2*numdof+2]=+kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[2*numdof+6]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[6*numdof+2]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[6*numdof+6]=+kmax*pow((IssmDouble)10.0,penalty_offset); + + Ke->values[3*numdof+3]=+kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[3*numdof+7]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[7*numdof+3]=-kmax*pow((IssmDouble)10.0,penalty_offset); + Ke->values[7*numdof+7]=+kmax*pow((IssmDouble)10.0,penalty_offset); /*Clean up and return*/ @@ -333,8 +333,8 @@ /*}}}*/ /*FUNCTION Penpair::PenaltyCreateKMatrixPrognostic {{{*/ -ElementMatrix* Penpair::PenaltyCreateKMatrixPrognostic(double kmax){ +ElementMatrix* Penpair::PenaltyCreateKMatrixPrognostic(IssmDouble kmax){ const int numdof=NUMVERTICES*NDOF1; - double penalty_factor; + IssmDouble penalty_factor; /*Initialize Element vector and return if necessary*/ @@ -345,8 +345,8 @@ //Create elementary matrix: add penalty to - Ke->values[0*numdof+0]=+kmax*pow((double)10.0,penalty_factor); - Ke->values[0*numdof+1]=-kmax*pow((double)10.0,penalty_factor); - Ke->values[1*numdof+0]=-kmax*pow((double)10.0,penalty_factor); - Ke->values[1*numdof+1]=+kmax*pow((double)10.0,penalty_factor); + Ke->values[0*numdof+0]=+kmax*pow((IssmDouble)10.0,penalty_factor); + Ke->values[0*numdof+1]=-kmax*pow((IssmDouble)10.0,penalty_factor); + Ke->values[1*numdof+0]=-kmax*pow((IssmDouble)10.0,penalty_factor); + Ke->values[1*numdof+1]=+kmax*pow((IssmDouble)10.0,penalty_factor); /*Clean up and return*/ Index: /issm/trunk-jpl/src/c/objects/Loads/Penpair.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Penpair.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Penpair.h (revision 12472) @@ -41,15 +41,15 @@ /*}}}*/ /*Update virtual functions resolution: {{{*/ - void InputUpdateFromVector(double* vector, int name, int type); + void InputUpdateFromVector(IssmDouble* vector, int name, int type); void InputUpdateFromVector(int* vector, int name, int type); void InputUpdateFromVector(bool* vector, int name, int type); - void InputUpdateFromMatrixDakota(double* matrix, int nrow, int ncols,int name, int type){_error_("Not implemented yet!");} - void InputUpdateFromVectorDakota(double* vector, int name, int type){_error_("Not implemented yet!");} + void InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrow, int ncols,int name, int type){_error_("Not implemented yet!");} + void InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){_error_("Not implemented yet!");} void InputUpdateFromVectorDakota(int* vector, int name, int type){_error_("Not implemented yet!");} void InputUpdateFromVectorDakota(bool* vector, int name, int type){_error_("Not implemented yet!");} - void InputUpdateFromConstant(double constant, int name); + void InputUpdateFromConstant(IssmDouble constant, int name); void InputUpdateFromConstant(int constant, int name); void InputUpdateFromConstant(bool constant, int name); - void InputUpdateFromSolution(double* solution){_error_("Not implemented yet!");} + void InputUpdateFromSolution(IssmDouble* solution){_error_("Not implemented yet!");} void InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");}; /*}}}*/ @@ -60,14 +60,14 @@ void CreatePVector(Vector* pf); void CreateJacobianMatrix(Matrix* Jff); - void PenaltyCreateKMatrix(Matrix* Kff,Matrix* Kfs,double kmax); - void PenaltyCreatePVector(Vector* pf, double kmax); - void PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax); + void PenaltyCreateKMatrix(Matrix* Kff,Matrix* Kfs,IssmDouble kmax); + void PenaltyCreatePVector(Vector* pf, IssmDouble kmax); + void PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax); bool InAnalysis(int analysis_type); /*}}}*/ /*Penpair management: {{{*/ - ElementMatrix* PenaltyCreateKMatrixDiagnosticHoriz(double kmax); - ElementMatrix* PenaltyCreateKMatrixDiagnosticMacAyealPattyn(double kmax); - ElementMatrix* PenaltyCreateKMatrixDiagnosticStokes(double kmax); - ElementMatrix* PenaltyCreateKMatrixPrognostic(double kmax); + ElementMatrix* PenaltyCreateKMatrixDiagnosticHoriz(IssmDouble kmax); + ElementMatrix* PenaltyCreateKMatrixDiagnosticMacAyealPattyn(IssmDouble kmax); + ElementMatrix* PenaltyCreateKMatrixDiagnosticStokes(IssmDouble kmax); + ElementMatrix* PenaltyCreateKMatrixPrognostic(IssmDouble kmax); /*}}}*/ }; Index: /issm/trunk-jpl/src/c/objects/Loads/Riftfront.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Riftfront.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Riftfront.cpp (revision 12472) @@ -45,6 +45,6 @@ int riftfront_type; int riftfront_fill; - double riftfront_friction; - double riftfront_fractionincrement; + IssmDouble riftfront_friction; + IssmDouble riftfront_fractionincrement; bool riftfront_shelf; int numberofelements; @@ -135,5 +135,5 @@ Input* input=NULL; int fill; - double friction,fractionincrement; + IssmDouble friction,fractionincrement; @@ -258,6 +258,6 @@ } /*}}}*/ -/*FUNCTION Riftfront::InputUpdateFromConstant(double constant,int name) {{{*/ -void Riftfront::InputUpdateFromConstant(double constant,int name){ +/*FUNCTION Riftfront::InputUpdateFromConstant(IssmDouble constant,int name) {{{*/ +void Riftfront::InputUpdateFromConstant(IssmDouble constant,int name){ /*Check that name is a Riftfront input*/ @@ -269,6 +269,6 @@ } /*}}}*/ -/*FUNCTION Riftfront::InputUpdateFromConstant(double* constant,int name) {{{*/ -void Riftfront::InputUpdateFromVector(double* vector, int name, int type){ +/*FUNCTION Riftfront::InputUpdateFromConstant(IssmDouble* constant,int name) {{{*/ +void Riftfront::InputUpdateFromVector(IssmDouble* vector, int name, int type){ /*Check that name is a Riftfront input*/ @@ -309,5 +309,5 @@ /*}}}*/ /*FUNCTION Riftfront::PenaltyCreateKMatrix {{{*/ -void Riftfront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){ +void Riftfront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,IssmDouble kmax){ /*Retrieve parameters: */ @@ -335,5 +335,5 @@ /*}}}*/ /*FUNCTION Riftfront::PenaltyCreatePVector {{{*/ -void Riftfront::PenaltyCreatePVector(Vector* pf,double kmax){ +void Riftfront::PenaltyCreatePVector(Vector* pf,IssmDouble kmax){ /*Retrieve parameters: */ @@ -381,14 +381,14 @@ /*Riftfront numerics*/ /*FUNCTION Riftfront::PenaltyCreateKMatrixDiagnosticHoriz {{{*/ -ElementMatrix* Riftfront::PenaltyCreateKMatrixDiagnosticHoriz(double kmax){ +ElementMatrix* Riftfront::PenaltyCreateKMatrixDiagnosticHoriz(IssmDouble kmax){ const int numdof = NDOF2*NUMVERTICES; int i,j; int dofs[1] = {0}; - double Ke_gg[4][4]; - double thickness; - double h[2]; - double penalty_offset; - double friction; + IssmDouble Ke_gg[4][4]; + IssmDouble thickness; + IssmDouble h[2]; + IssmDouble penalty_offset; + IssmDouble friction; /*Objects: */ @@ -464,20 +464,20 @@ /*}}}*/ /*FUNCTION Riftfront::PenaltyCreatePVectorDiagnosticHoriz {{{*/ -ElementVector* Riftfront::PenaltyCreatePVectorDiagnosticHoriz(double kmax){ +ElementVector* Riftfront::PenaltyCreatePVectorDiagnosticHoriz(IssmDouble kmax){ const int numdof = NDOF2*NUMVERTICES; int i,j; - double rho_ice; - double rho_water; - double gravity; - double thickness; - double h[2]; - double bed; - double b[2]; - double pressure; - double pressure_litho; - double pressure_air; - double pressure_melange; - double pressure_water; + IssmDouble rho_ice; + IssmDouble rho_water; + IssmDouble gravity; + IssmDouble thickness; + IssmDouble h[2]; + IssmDouble bed; + IssmDouble b[2]; + IssmDouble pressure; + IssmDouble pressure_litho; + IssmDouble pressure_air; + IssmDouble pressure_melange; + IssmDouble pressure_water; int fill; bool shelf; @@ -521,13 +521,13 @@ if(shelf){ /*We are on an ice shelf, hydrostatic equilibrium is used to determine the pressure for water fill: */ - pressure=rho_ice*gravity*pow(thickness,(double)2)/(double)2 - rho_water*gravity*pow(bed,(double)2)/(double)2; + pressure=rho_ice*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2 - rho_water*gravity*pow(bed,(IssmDouble)2)/(IssmDouble)2; } else{ //We are on an icesheet, we assume the water column fills the entire front: */ - pressure=rho_ice*gravity*pow(thickness,(double)2)/(double)2 - rho_water*gravity*pow(thickness,(double)2)/(double)2; + pressure=rho_ice*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2 - rho_water*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2; } } else if(fill==AirEnum){ - pressure=rho_ice*gravity*pow(thickness,(double)2)/(double)2; //icefront on an ice sheet, pressure imbalance ice vs air. + pressure=rho_ice*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2; //icefront on an ice sheet, pressure imbalance ice vs air. } else if(fill==IceEnum){ //icefront finding itself against another icefront (pressure imbalance is fully compensated, ice vs ice) @@ -538,7 +538,7 @@ if(!shelf) _error_("%s%i%s","fill type ",fill," not supported on ice sheets yet."); - pressure_litho=rho_ice*gravity*pow(thickness,(double)2)/(double)2; + pressure_litho=rho_ice*gravity*pow(thickness,(IssmDouble)2)/(IssmDouble)2; pressure_air=0; - pressure_melange=rho_ice*gravity*pow(fraction*thickness,(double)2)/(double)2; + pressure_melange=rho_ice*gravity*pow(fraction*thickness,(IssmDouble)2)/(IssmDouble)2; pressure_water=1.0/2.0*rho_water*gravity* ( pow(bed,2.0)-pow(rho_ice/rho_water*fraction*thickness,2.0) ); @@ -569,14 +569,14 @@ const int numnodes = 2; - double max_penetration; - double penetration; + IssmDouble max_penetration; + IssmDouble penetration; int activate; int found; int unstable; - double vx1; - double vy1; - double vx2; - double vy2; - double fractionincrement; + IssmDouble vx1; + IssmDouble vy1; + IssmDouble vx2; + IssmDouble vy2; + IssmDouble fractionincrement; /*Objects: */ @@ -632,5 +632,5 @@ /*increase melange fraction: */ this->fraction+=fractionincrement; - if (this->fraction>1)this->fraction=(double)1.0; + if (this->fraction>1)this->fraction=(IssmDouble)1.0; //printf("riftfront %i fraction: %g\n",this->Id(),this->fraction); } @@ -675,5 +675,5 @@ int found=0; - double converged=0; + IssmDouble converged=0; this->inputs->GetInputValue(&converged,ConvergedEnum); @@ -690,14 +690,14 @@ /*}}}*/ /*FUNCTION Riftfront::MaxPenetration {{{*/ -int Riftfront::MaxPenetration(double* ppenetration){ +int Riftfront::MaxPenetration(IssmDouble* ppenetration){ const int numnodes=2; - double max_penetration; - double penetration=0; + IssmDouble max_penetration; + IssmDouble penetration=0; int found; - double vx1; - double vy1; - double vx2; - double vy2; + IssmDouble vx1; + IssmDouble vy1; + IssmDouble vx2; + IssmDouble vy2; /*Objects: */ @@ -736,12 +736,12 @@ /*}}}*/ /*FUNCTION Riftfront::Penetration {{{*/ -int Riftfront::Penetration(double* ppenetration){ - - double vx1; - double vy1; - double vx2; - double vy2; - - double penetration; +int Riftfront::Penetration(IssmDouble* ppenetration){ + + IssmDouble vx1; + IssmDouble vy1; + IssmDouble vx2; + IssmDouble vy2; + + IssmDouble penetration; int found; @@ -779,13 +779,13 @@ const int numnodes = 2; - double max_penetration; - double penetration; + IssmDouble max_penetration; + IssmDouble penetration; int activate; int unstable; int found; - double vx1; - double vy1; - double vx2; - double vy2; + IssmDouble vx1; + IssmDouble vy1; + IssmDouble vx2; + IssmDouble vy2; /*Objects: */ @@ -831,11 +831,11 @@ const int numnodes = 2; - double penetration; + IssmDouble penetration; int unstable; int found; - double vx1; - double vy1; - double vx2; - double vy2; + IssmDouble vx1; + IssmDouble vy1; + IssmDouble vx2; + IssmDouble vy2; /*Objects: */ Index: /issm/trunk-jpl/src/c/objects/Loads/Riftfront.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Loads/Riftfront.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Loads/Riftfront.h (revision 12472) @@ -38,7 +38,7 @@ bool prestable; bool material_converged; - double normal[2]; - double length; - double fraction; + IssmDouble normal[2]; + IssmDouble length; + IssmDouble fraction; int state; @@ -61,15 +61,15 @@ /*}}}*/ /*Update virtual functions resolution: {{{*/ - void InputUpdateFromVector(double* vector, int name, int type); + void InputUpdateFromVector(IssmDouble* vector, int name, int type); void InputUpdateFromVector(int* vector, int name, int type){_error_("Not implemented yet!");} void InputUpdateFromVector(bool* vector, int name, int type){_error_("Not implemented yet!");} - void InputUpdateFromMatrixDakota(double* matrix, int nrows,int ncols, int name, int type){_error_("Not implemented yet!");} - void InputUpdateFromVectorDakota(double* vector, int name, int type){_error_("Not implemented yet!");} + void InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows,int ncols, int name, int type){_error_("Not implemented yet!");} + void InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){_error_("Not implemented yet!");} void InputUpdateFromVectorDakota(int* vector, int name, int type){_error_("Not implemented yet!");} void InputUpdateFromVectorDakota(bool* vector, int name, int type){_error_("Not implemented yet!");} - void InputUpdateFromConstant(double constant, int name); + void InputUpdateFromConstant(IssmDouble constant, int name); void InputUpdateFromConstant(int constant, int name){_error_("Not implemented yet!");} void InputUpdateFromConstant(bool constant, int name); - void InputUpdateFromSolution(double* solution){_error_("Not implemented yet!");} + void InputUpdateFromSolution(IssmDouble* solution){_error_("Not implemented yet!");} void InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");}; /*}}}*/ @@ -80,13 +80,13 @@ void CreatePVector(Vector* pf); void CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");}; - void PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");}; - void PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax); - void PenaltyCreatePVector(Vector* pf, double kmax); + void PenaltyCreateJacobianMatrix(Matrix* Jff,IssmDouble kmax){_error_("Not implemented yet");}; + void PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, IssmDouble kmax); + void PenaltyCreatePVector(Vector* pf, IssmDouble kmax); bool InAnalysis(int analysis_type); /*}}}*/ /*Riftfront specific routines: {{{*/ bool PreStable(); - ElementMatrix* PenaltyCreateKMatrixDiagnosticHoriz(double kmax); - ElementVector* PenaltyCreatePVectorDiagnosticHoriz(double kmax); + ElementMatrix* PenaltyCreateKMatrixDiagnosticHoriz(IssmDouble kmax); + ElementVector* PenaltyCreatePVectorDiagnosticHoriz(IssmDouble kmax); void SetPreStable(); int PreConstrain(int* punstable); @@ -94,6 +94,6 @@ void FreezeConstraints(void); bool IsFrozen(void); - int Penetration(double* ppenetration); - int MaxPenetration(double* ppenetration); + int Penetration(IssmDouble* ppenetration); + int MaxPenetration(IssmDouble* ppenetration); int PotentialUnstableConstraint(int* punstable); int IsMaterialStable(void); Index: /issm/trunk-jpl/src/c/objects/Materials/Matice.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Materials/Matice.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Materials/Matice.cpp (revision 12472) @@ -129,8 +129,8 @@ /*}}}*/ /*FUNCTION Matice::GetB {{{*/ -double Matice::GetB(){ +IssmDouble Matice::GetB(){ /*Output*/ - double B; + IssmDouble B; inputs->GetInputAverage(&B,MaterialsRheologyBEnum); @@ -139,8 +139,8 @@ /*}}}*/ /*FUNCTION Matice::GetBbar {{{*/ -double Matice::GetBbar(){ +IssmDouble Matice::GetBbar(){ /*Output*/ - double Bbar; + IssmDouble Bbar; inputs->GetInputAverage(&Bbar,MaterialsRheologyBbarEnum); @@ -149,8 +149,8 @@ /*}}}*/ /*FUNCTION Matice::GetN {{{*/ -double Matice::GetN(){ +IssmDouble Matice::GetN(){ /*Output*/ - double n; + IssmDouble n; inputs->GetInputAverage(&n,MaterialsRheologyNEnum); @@ -191,5 +191,5 @@ /*}}}*/ /*FUNCTION Matice::GetViscosity2d {{{*/ -void Matice::GetViscosity2d(double* pviscosity, double* epsilon){ +void Matice::GetViscosity2d(IssmDouble* pviscosity, IssmDouble* epsilon){ /*From a string tensor and a material object, return viscosity, using Glen's flow law. B @@ -205,12 +205,12 @@ /*output: */ - double viscosity; + IssmDouble viscosity; /*input strain rate: */ - double exx,eyy,exy; + IssmDouble exx,eyy,exy; /*Intermediary: */ - double A,e; - double B,n; + IssmDouble A,e; + IssmDouble B,n; /*Get B and n*/ @@ -224,5 +224,5 @@ else{ if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0)){ - viscosity=0.5*pow((double)10,(double)14); + viscosity=0.5*pow((IssmDouble)10,(IssmDouble)14); } else{ @@ -236,5 +236,5 @@ if(A==0){ /*Maxiviscositym viscosity for 0 shear areas: */ - viscosity=2.5*pow((double)10,(double)17); + viscosity=2.5*pow((IssmDouble)10,(IssmDouble)17); } else{ @@ -255,5 +255,5 @@ /*}}}*/ /*FUNCTION Matice::GetViscosity3d {{{*/ -void Matice::GetViscosity3d(double* pviscosity3d, double* epsilon){ +void Matice::GetViscosity3d(IssmDouble* pviscosity3d, IssmDouble* epsilon){ /*Return viscosity accounting for steady state power law creep [Thomas and MacAyeal, 1982]: @@ -271,12 +271,12 @@ /*output: */ - double viscosity3d; + IssmDouble viscosity3d; /*input strain rate: */ - double exx,eyy,exy,exz,eyz; + IssmDouble exx,eyy,exy,exz,eyz; /*Intermediaries: */ - double A,e; - double B,n; + IssmDouble A,e; + IssmDouble B,n; /*Get B and n*/ @@ -291,5 +291,5 @@ if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0) && (epsilon[3]==0) && (epsilon[4]==0)){ - viscosity3d=0.5*pow((double)10,(double)14); + viscosity3d=0.5*pow((IssmDouble)10,(IssmDouble)14); } else{ @@ -306,5 +306,5 @@ if(A==0){ /*Maxiviscosity3dm viscosity for 0 shear areas: */ - viscosity3d=2.25*pow((double)10,(double)17); + viscosity3d=2.25*pow((IssmDouble)10,(IssmDouble)17); } else{ @@ -326,5 +326,5 @@ /*}}}*/ /*FUNCTION Matice::GetViscosity3dStokes {{{*/ -void Matice::GetViscosity3dStokes(double* pviscosity3d, double* epsilon){ +void Matice::GetViscosity3dStokes(IssmDouble* pviscosity3d, IssmDouble* epsilon){ /*Return viscosity accounting for steady state power law creep [Thomas and MacAyeal, 1982]: * @@ -341,16 +341,16 @@ /*output: */ - double viscosity3d; + IssmDouble viscosity3d; /*input strain rate: */ - double exx,eyy,exy,exz,eyz,ezz; + IssmDouble exx,eyy,exy,exz,eyz,ezz; /*Intermediaries: */ - double A,e; - double B,n; - double eps0; + IssmDouble A,e; + IssmDouble B,n; + IssmDouble eps0; /*Get B and n*/ - eps0=pow((double)10,(double)-27); + eps0=pow((IssmDouble)10,(IssmDouble)-27); B=GetB(); n=GetN(); @@ -363,5 +363,5 @@ if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0) && (epsilon[3]==0) && (epsilon[4]==0) && (epsilon[5]==0)){ - viscosity3d=0.5*pow((double)10,(double)14); + viscosity3d=0.5*pow((IssmDouble)10,(IssmDouble)14); } else{ @@ -379,5 +379,5 @@ if(A==0){ /*Maxiviscosity3dm viscosity for 0 shear areas: */ - viscosity3d=2.25*pow((double)10,(double)17); + viscosity3d=2.25*pow((IssmDouble)10,(IssmDouble)17); } else{ @@ -398,5 +398,5 @@ /*}}}*/ /*FUNCTION Matice::GetViscosityComplement {{{*/ -void Matice::GetViscosityComplement(double* pviscosity_complement, double* epsilon){ +void Matice::GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* epsilon){ /*Return viscosity accounting for steady state power law creep [Thomas and MacAyeal, 1982]: * @@ -410,12 +410,12 @@ /*output: */ - double viscosity_complement; + IssmDouble viscosity_complement; /*input strain rate: */ - double exx,eyy,exy; + IssmDouble exx,eyy,exy; /*Intermediary value A and exponent e: */ - double A,e; - double B,n; + IssmDouble A,e; + IssmDouble B,n; /*Get B and n*/ @@ -432,5 +432,5 @@ if(A==0){ /*Maximum viscosity_complement for 0 shear areas: */ - viscosity_complement=2.25*pow((double)10,(double)17); + viscosity_complement=2.25*pow((IssmDouble)10,(IssmDouble)17); } else{ @@ -441,5 +441,5 @@ } else{ - viscosity_complement=4.5*pow((double)10,(double)17); + viscosity_complement=4.5*pow((IssmDouble)10,(IssmDouble)17); } @@ -454,12 +454,12 @@ /*}}}*/ /*FUNCTION Matice::GetViscosityDerivativeEpsSquare{{{*/ -void Matice::GetViscosityDerivativeEpsSquare(double* pmu_prime, double* epsilon){ +void Matice::GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* epsilon){ /*output: */ - double mu_prime; - double mu,n,eff2; + IssmDouble mu_prime; + IssmDouble mu,n,eff2; /*input strain rate: */ - double exx,eyy,exy,exz,eyz; + IssmDouble exx,eyy,exy,exz,eyz; /*Get visocisty and n*/ @@ -469,5 +469,5 @@ if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0) && (epsilon[3]==0) && (epsilon[4]==0)){ - mu_prime=0.5*pow((double)10,(double)14); + mu_prime=0.5*pow((IssmDouble)10,(IssmDouble)14); } else{ @@ -488,12 +488,12 @@ /*}}}*/ /*FUNCTION Matice::GetViscosity2dDerivativeEpsSquare{{{*/ -void Matice::GetViscosity2dDerivativeEpsSquare(double* pmu_prime, double* epsilon){ +void Matice::GetViscosity2dDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* epsilon){ /*output: */ - double mu_prime; - double mu,n,eff2; + IssmDouble mu_prime; + IssmDouble mu,n,eff2; /*input strain rate: */ - double exx,eyy,exy,exz; + IssmDouble exx,eyy,exy,exz; /*Get visocisty and n*/ @@ -502,5 +502,5 @@ if((epsilon[0]==0) && (epsilon[1]==0) && (epsilon[2]==0)){ - mu_prime=0.5*pow((double)10,(double)14); + mu_prime=0.5*pow((IssmDouble)10,(IssmDouble)14); } else{ @@ -526,6 +526,6 @@ } /*}}}*/ -/*FUNCTION Matice::InputUpdateFromVector(double* vector, int name, int type) {{{*/ -void Matice::InputUpdateFromVector(double* vector, int name, int type){ +/*FUNCTION Matice::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/ +void Matice::InputUpdateFromVector(IssmDouble* vector, int name, int type){ /*Intermediaries*/ @@ -545,5 +545,5 @@ case TriaEnum: - double values[3]; + IssmDouble values[3]; for (int i=0;i<3;i++) values[i]=vector[((Tria*)element)->nodes[i]->GetVertexDof()]; this->inputs->AddInput(new TriaP1Input(name,values)); @@ -566,6 +566,6 @@ } /*}}}*/ -/*FUNCTION Matice::InputUpdateFromVectorDakota(double* vector, int name, int type) {{{*/ -void Matice::InputUpdateFromVectorDakota(double* vector, int name, int type){ +/*FUNCTION Matice::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type) {{{*/ +void Matice::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){ /*Intermediaries*/ @@ -587,5 +587,5 @@ case TriaEnum: - double values[3]; + IssmDouble values[3]; for (int i=0;i<3;i++) values[i]=vector[((Tria*)element)->nodes[i]->GetSidList()]; //use sid list, to index into serial oriented vector this->inputs->AddInput(new TriaP1Input(name,values)); @@ -618,5 +618,5 @@ /*}}}*/ /*FUNCTION Matice::InputUpdateFromMatrixDakota(int* vector, int name, int type) {{{*/ -void Matice::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols,int name, int type){ +void Matice::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols,int name, int type){ /*Nothing updated yet*/ } @@ -632,6 +632,6 @@ } /*}}}*/ -/*FUNCTION Matice::InputUpdateFromConstant(double constant, int name) {{{*/ -void Matice::InputUpdateFromConstant(double constant, int name){ +/*FUNCTION Matice::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/ +void Matice::InputUpdateFromConstant(IssmDouble constant, int name){ /*Nothing updated yet*/ } @@ -648,5 +648,5 @@ /*}}}*/ /*FUNCTION Matice::InputUpdateFromSolution{{{*/ -void Matice::InputUpdateFromSolution(double* solution){ +void Matice::InputUpdateFromSolution(IssmDouble* solution){ /*Nothing updated yet*/ } @@ -671,7 +671,7 @@ /*Intermediaries*/ const int num_vertices = 3; //Tria has 3 vertices - double nodeinputs[num_vertices]; - double cmmininputs[num_vertices]; - double cmmaxinputs[num_vertices]; + IssmDouble nodeinputs[num_vertices]; + IssmDouble cmmininputs[num_vertices]; + IssmDouble cmmaxinputs[num_vertices]; /*Get B*/ @@ -713,7 +713,7 @@ /*Intermediaries*/ const int num_vertices = 6; //Penta has 6 vertices - double nodeinputs[num_vertices]; - double cmmininputs[num_vertices]; - double cmmaxinputs[num_vertices]; + IssmDouble nodeinputs[num_vertices]; + IssmDouble cmmininputs[num_vertices]; + IssmDouble cmmaxinputs[num_vertices]; /*Get B*/ Index: /issm/trunk-jpl/src/c/objects/Materials/Matice.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Materials/Matice.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Materials/Matice.h (revision 12472) @@ -39,15 +39,15 @@ /*}}}*/ /*Update virtual functions definitions: {{{*/ - void InputUpdateFromVector(double* vector, int name, int type); + void InputUpdateFromVector(IssmDouble* vector, int name, int type); void InputUpdateFromVector(int* vector, int name, int type); void InputUpdateFromVector(bool* vector, int name, int type); - void InputUpdateFromMatrixDakota(double* matrix, int nrow, int ncols, int name, int type); - void InputUpdateFromVectorDakota(double* vector, int name, int type); + void InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrow, int ncols, int name, int type); + void InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type); void InputUpdateFromVectorDakota(int* vector, int name, int type); void InputUpdateFromVectorDakota(bool* vector, int name, int type); - void InputUpdateFromConstant(double constant, int name); + void InputUpdateFromConstant(IssmDouble constant, int name); void InputUpdateFromConstant(int constant, int name); void InputUpdateFromConstant(bool constant, int name); - void InputUpdateFromSolution(double* solution); + void InputUpdateFromSolution(IssmDouble* solution); void InputUpdateFromIoModel(int index, IoModel* iomodel); /*}}}*/ @@ -59,13 +59,13 @@ /*Matice Numerics: {{{*/ void SetCurrentConfiguration(Elements* elementsin,Loads* loadsin,Nodes* nodesin,Vertices* verticesin,Materials* materialsin,Parameters* parametersin); - void GetViscosity2d(double* pviscosity, double* pepsilon); - void GetViscosity3d(double* pviscosity3d, double* pepsilon); - void GetViscosity3dStokes(double* pviscosity3d, double* epsilon); - void GetViscosityComplement(double* pviscosity_complement, double* pepsilon); - void GetViscosityDerivativeEpsSquare(double* pmu_prime, double* pepsilon); - void GetViscosity2dDerivativeEpsSquare(double* pmu_prime, double* pepsilon); - double GetB(); - double GetBbar(); - double GetN(); + void GetViscosity2d(IssmDouble* pviscosity, IssmDouble* pepsilon); + void GetViscosity3d(IssmDouble* pviscosity3d, IssmDouble* pepsilon); + void GetViscosity3dStokes(IssmDouble* pviscosity3d, IssmDouble* epsilon); + void GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon); + void GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon); + void GetViscosity2dDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon); + IssmDouble GetB(); + IssmDouble GetBbar(); + IssmDouble GetN(); bool IsInput(int name); /*}}}*/ Index: /issm/trunk-jpl/src/c/objects/Materials/Matpar.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Materials/Matpar.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Materials/Matpar.cpp (revision 12472) @@ -104,6 +104,6 @@ /*Update virtual functions definitions:*/ -/*FUNCTION Matpar::InputUpdateFromVector(double* vector, int name, int type) {{{*/ -void Matpar::InputUpdateFromVector(double* vector, int name, int type){ +/*FUNCTION Matpar::InputUpdateFromVector(IssmDouble* vector, int name, int type) {{{*/ +void Matpar::InputUpdateFromVector(IssmDouble* vector, int name, int type){ /*Nothing updated yet*/ } @@ -119,6 +119,6 @@ } /*}}}*/ -/*FUNCTION Matpar::InputUpdateFromVectorDakota(double* vector, int name, int type) {{{*/ -void Matpar::InputUpdateFromVectorDakota(double* vector, int name, int type){ +/*FUNCTION Matpar::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type) {{{*/ +void Matpar::InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){ /*Nothing updated yet*/ } @@ -135,10 +135,10 @@ /*}}}*/ /*FUNCTION Matpar::InputUpdateFromMatrixDakota(int* vector, int name, int type) {{{*/ -void Matpar::InputUpdateFromMatrixDakota(double* matrix, int nrows, int ncols,int name, int type){ - /*Nothing updated yet*/ -} -/*}}}*/ -/*FUNCTION Matpar::InputUpdateFromConstant(double constant, int name) {{{*/ -void Matpar::InputUpdateFromConstant(double constant, int name){ +void Matpar::InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols,int name, int type){ + /*Nothing updated yet*/ +} +/*}}}*/ +/*FUNCTION Matpar::InputUpdateFromConstant(IssmDouble constant, int name) {{{*/ +void Matpar::InputUpdateFromConstant(IssmDouble constant, int name){ switch(name){ @@ -199,5 +199,5 @@ /*}}}*/ /*FUNCTION Matpar::InputUpdateFromSolution{{{*/ -void Matpar::InputUpdateFromSolution(double* solution){ +void Matpar::InputUpdateFromSolution(IssmDouble* solution){ /*Nothing updated yet*/ } @@ -213,40 +213,40 @@ /*}}}*/ /*FUNCTION Matpar::GetBeta {{{*/ -double Matpar::GetBeta(){ +IssmDouble Matpar::GetBeta(){ return beta; } /*}}}*/ /*FUNCTION Matpar::GetG {{{*/ -double Matpar::GetG(){ +IssmDouble Matpar::GetG(){ return g; } /*}}}*/ /*FUNCTION Matpar::GetHeatCapacity {{{*/ -double Matpar::GetHeatCapacity(){ +IssmDouble Matpar::GetHeatCapacity(){ return heatcapacity; } /*}}}*/ /*FUNCTION Matpar::GetLatentHeat {{{*/ -double Matpar::GetLatentHeat(){ +IssmDouble Matpar::GetLatentHeat(){ return latentheat; } /*}}}*/ /*FUNCTION Matpar::GetMeltingPoint {{{*/ -double Matpar::GetMeltingPoint(){ +IssmDouble Matpar::GetMeltingPoint(){ return meltingpoint; } /*}}}*/ /*FUNCTION Matpar::GetReferenceTemperature {{{*/ -double Matpar::GetReferenceTemperature(){ +IssmDouble Matpar::GetReferenceTemperature(){ return referencetemperature; } /*}}}*/ /*FUNCTION Matpar::GetMixedLayerCapacity {{{*/ -double Matpar::GetMixedLayerCapacity(){ +IssmDouble Matpar::GetMixedLayerCapacity(){ return mixed_layer_capacity; } /*}}}*/ /*FUNCTION Matpar::GetRhoIce {{{*/ -double Matpar::GetRhoIce(){ +IssmDouble Matpar::GetRhoIce(){ return rho_ice; @@ -254,65 +254,65 @@ /*}}}*/ /*FUNCTION Matpar::GetRhoWater {{{*/ -double Matpar::GetRhoWater(){ +IssmDouble Matpar::GetRhoWater(){ return rho_water; } /*}}}*/ /*FUNCTION Matpar::GetRhoFreshwater {{{*/ -double Matpar::GetRhoFreshwater(){ +IssmDouble Matpar::GetRhoFreshwater(){ return rho_freshwater; } /*}}}*/ /*FUNCTION Matpar::GetMuWater {{{*/ -double Matpar::GetMuWater(){ +IssmDouble Matpar::GetMuWater(){ return mu_water; } /*}}}*/ /*FUNCTION Matpar::GetThermalConductivity {{{*/ -double Matpar::GetThermalConductivity(){ +IssmDouble Matpar::GetThermalConductivity(){ return thermalconductivity; } /*}}}*/ /*FUNCTION Matpar::GetThermalExchangeVelocity {{{*/ -double Matpar::GetThermalExchangeVelocity(){ +IssmDouble Matpar::GetThermalExchangeVelocity(){ return thermal_exchange_velocity; } /*}}}*/ /*FUNCTION Matpar::GetKn {{{*/ -double Matpar::GetKn(){ +IssmDouble Matpar::GetKn(){ return kn; } /*}}}*/ /*FUNCTION Matpar::GetHydrologyP {{{*/ -double Matpar::GetHydrologyP(){ +IssmDouble Matpar::GetHydrologyP(){ return hydro_p; } /*}}}*/ /*FUNCTION Matqar::GetHydrologyQ {{{*/ -double Matpar::GetHydrologyQ(){ +IssmDouble Matpar::GetHydrologyQ(){ return hydro_q; } /*}}}*/ /*FUNCTION Matpar::GetHydrologyCR {{{*/ -double Matpar::GetHydrologyCR(){ +IssmDouble Matpar::GetHydrologyCR(){ return hydro_CR; } /*}}}*/ /*FUNCTION Matpar::GetHydrologyN {{{*/ -double Matpar::GetHydrologyN(){ +IssmDouble Matpar::GetHydrologyN(){ return hydro_n; } /*}}}*/ /*FUNCTION Matpar::TMeltingPoint {{{*/ -double Matpar::TMeltingPoint(double pressure){ +IssmDouble Matpar::TMeltingPoint(IssmDouble pressure){ return meltingpoint-beta*pressure; } /*}}}*/ /*FUNCTION Matpar::PureIceEnthalpy{{{*/ -double Matpar::PureIceEnthalpy(double pressure){ +IssmDouble Matpar::PureIceEnthalpy(IssmDouble pressure){ return heatcapacity*(TMeltingPoint(pressure)-referencetemperature); } /*}}}*/ /*FUNCTION Matpar::GetEnthalpyDiffusionParameter{{{*/ -double Matpar::GetEnthalpyDiffusionParameter(double enthalpy,double pressure){ +IssmDouble Matpar::GetEnthalpyDiffusionParameter(IssmDouble enthalpy,IssmDouble pressure){ if(enthalpygglobaldoflist=xNew(this->nrows); - this->values=xNewZeroInit(this->nrows*this->ncols); + this->values=xNewZeroInit(this->nrows*this->ncols); for(i=0;inrows;i++){ for(j=0;jncols;j++){ @@ -205,5 +205,5 @@ /*fill values with 0: */ - this->values=xNewZeroInit(this->nrows*this->ncols); + this->values=xNewZeroInit(this->nrows*this->ncols); /*g list*/ @@ -231,5 +231,5 @@ ElementMatrix::~ElementMatrix(){ - xDelete(this->values); + xDelete(this->values); xDelete(this->gglobaldoflist); xDelete(this->row_flocaldoflist); @@ -249,5 +249,5 @@ int i,j; - double* localvalues=NULL; + IssmDouble* localvalues=NULL; /*If Kfs is not provided, call the other function*/ @@ -265,5 +265,5 @@ if(this->row_fsize){ /*first, retrieve values that are in the f-set from the g-set values matrix: */ - localvalues=xNew(this->row_fsize*this->row_fsize); + localvalues=xNew(this->row_fsize*this->row_fsize); for(i=0;irow_fsize;i++){ for(j=0;jrow_fsize;j++){ @@ -275,5 +275,5 @@ /*Free ressources:*/ - xDelete(localvalues); + xDelete(localvalues); } @@ -281,5 +281,5 @@ if((this->row_ssize!=0) && (this->row_fsize!=0)){ /*first, retrieve values that are in the f and s-set from the g-set values matrix: */ - localvalues=xNew(this->row_fsize*this->row_ssize); + localvalues=xNew(this->row_fsize*this->row_ssize); for(i=0;irow_fsize;i++){ for(j=0;jrow_ssize;j++){ @@ -291,5 +291,5 @@ /*Free ressources:*/ - xDelete(localvalues); + xDelete(localvalues); } } @@ -304,5 +304,5 @@ int i,j; - double* localvalues=NULL; + IssmDouble* localvalues=NULL; /*Check that Jff is not NULL*/ @@ -317,5 +317,5 @@ if(this->row_fsize){ /*first, retrieve values that are in the f-set from the g-set values matrix: */ - localvalues=xNew(this->row_fsize*this->row_fsize); + localvalues=xNew(this->row_fsize*this->row_fsize); for(i=0;irow_fsize;i++){ for(j=0;jrow_fsize;j++){ @@ -327,5 +327,5 @@ /*Free ressources:*/ - xDelete(localvalues); + xDelete(localvalues); } @@ -428,16 +428,16 @@ this->dofsymmetrical=Ke->dofsymmetrical; - this->values=xNew(this->nrows*this->ncols); - memcpy(this->values,Ke->values,this->nrows*this->ncols*sizeof(double)); + this->values=xNew(this->nrows*this->ncols); + xMemCpy(this->values,Ke->values,this->nrows*this->ncols); this->gglobaldoflist=xNew(this->nrows); - memcpy(this->gglobaldoflist,Ke->gglobaldoflist,this->nrows*sizeof(int)); + xMemCpy(this->gglobaldoflist,Ke->gglobaldoflist,this->nrows); this->row_fsize=Ke->row_fsize; if(this->row_fsize){ this->row_flocaldoflist=xNew(this->row_fsize); - memcpy(this->row_flocaldoflist,Ke->row_flocaldoflist,this->row_fsize*sizeof(int)); + xMemCpy(this->row_flocaldoflist,Ke->row_flocaldoflist,this->row_fsize); this->row_fglobaldoflist=xNew(this->row_fsize); - memcpy(this->row_fglobaldoflist,Ke->row_fglobaldoflist,this->row_fsize*sizeof(int)); + xMemCpy(this->row_fglobaldoflist,Ke->row_fglobaldoflist,this->row_fsize); } else{ @@ -449,7 +449,7 @@ if(this->row_ssize){ this->row_slocaldoflist=xNew(this->row_ssize); - memcpy(this->row_slocaldoflist,Ke->row_slocaldoflist,this->row_ssize*sizeof(int)); + xMemCpy(this->row_slocaldoflist,Ke->row_slocaldoflist,this->row_ssize); this->row_sglobaldoflist=xNew(this->row_ssize); - memcpy(this->row_sglobaldoflist,Ke->row_sglobaldoflist,this->row_ssize*sizeof(int)); + xMemCpy(this->row_sglobaldoflist,Ke->row_sglobaldoflist,this->row_ssize); } else{ @@ -461,7 +461,7 @@ if(this->col_fsize){ this->col_flocaldoflist=xNew(this->col_fsize); - memcpy(this->col_flocaldoflist,Ke->col_flocaldoflist,this->col_fsize*sizeof(int)); + xMemCpy(this->col_flocaldoflist,Ke->col_flocaldoflist,this->col_fsize); this->col_fglobaldoflist=xNew(this->col_fsize); - memcpy(this->col_fglobaldoflist,Ke->col_fglobaldoflist,this->col_fsize*sizeof(int)); + xMemCpy(this->col_fglobaldoflist,Ke->col_fglobaldoflist,this->col_fsize); } else{ @@ -473,7 +473,7 @@ if(this->col_ssize){ this->col_slocaldoflist=xNew(this->col_ssize); - memcpy(this->col_slocaldoflist,Ke->col_slocaldoflist,this->col_ssize*sizeof(int)); + xMemCpy(this->col_slocaldoflist,Ke->col_slocaldoflist,this->col_ssize); this->col_sglobaldoflist=xNew(this->col_ssize); - memcpy(this->col_sglobaldoflist,Ke->col_sglobaldoflist,this->col_ssize*sizeof(int)); + xMemCpy(this->col_sglobaldoflist,Ke->col_sglobaldoflist,this->col_ssize); } else{ @@ -484,5 +484,5 @@ /*}}}*/ /*FUNCTION ElementMatrix::SetDiag{{{*/ -void ElementMatrix::SetDiag(double scalar){ +void ElementMatrix::SetDiag(IssmDouble scalar){ int i; Index: /issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.h (revision 12472) @@ -24,5 +24,5 @@ int ncols; bool dofsymmetrical; - double* values; + IssmDouble* values; //gset @@ -64,5 +64,5 @@ void Transpose(void); void Init(ElementMatrix* Ke); - void SetDiag(double scalar); + void SetDiag(IssmDouble scalar); /*}}}*/ }; Index: /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp (revision 12472) @@ -75,5 +75,5 @@ /*Gset and values*/ this->gglobaldoflist=xNew(this->nrows); - this->values=xNewZeroInit(this->nrows); + this->values=xNewZeroInit(this->nrows); for(i=0;inrows;i++){ this->values[i] += pe1->values[i]; @@ -138,5 +138,5 @@ /*fill values with 0: */ - this->values=xNewZeroInit(this->nrows); + this->values=xNewZeroInit(this->nrows); /*g list*/ @@ -152,5 +152,5 @@ ElementVector::~ElementVector(){ - xDelete(this->values); + xDelete(this->values); xDelete(this->gglobaldoflist); xDelete(this->flocaldoflist); @@ -164,5 +164,5 @@ int i; - double* localvalues=NULL; + IssmDouble* localvalues=NULL; /*In debugging mode, check consistency (no NaN, and values not too big)*/ @@ -171,5 +171,5 @@ if(this->fsize){ /*first, retrieve values that are in the f-set from the g-set values vector: */ - localvalues=xNew(this->fsize); + localvalues=xNew(this->fsize); for(i=0;ifsize;i++){ localvalues[i]=this->values[this->flocaldoflist[i]]; @@ -179,5 +179,5 @@ /*Free ressources:*/ - xDelete(localvalues); + xDelete(localvalues); } @@ -188,9 +188,9 @@ int i; - double* localvalues=NULL; + IssmDouble* localvalues=NULL; if(this->fsize){ /*first, retrieve values that are in the f-set from the g-set values vector: */ - localvalues=xNew(this->fsize); + localvalues=xNew(this->fsize); for(i=0;ifsize;i++){ localvalues[i]=this->values[this->flocaldoflist[i]]; @@ -200,5 +200,5 @@ /*Free ressources:*/ - xDelete(localvalues); + xDelete(localvalues); } @@ -245,16 +245,16 @@ this->nrows =pe->nrows; - this->values=xNew(this->nrows); - memcpy(this->values,pe->values,this->nrows*sizeof(double)); + this->values=xNew(this->nrows); + xMemCpy(this->values,pe->values,this->nrows); this->gglobaldoflist=xNew(this->nrows); - memcpy(this->gglobaldoflist,pe->gglobaldoflist,this->nrows*sizeof(int)); + xMemCpy(this->gglobaldoflist,pe->gglobaldoflist,this->nrows); this->fsize=pe->fsize; if(this->fsize){ this->flocaldoflist=xNew(this->fsize); - memcpy(this->flocaldoflist,pe->flocaldoflist,this->fsize*sizeof(int)); + xMemCpy(this->flocaldoflist,pe->flocaldoflist,this->fsize); this->fglobaldoflist=xNew(this->fsize); - memcpy(this->fglobaldoflist,pe->fglobaldoflist,this->fsize*sizeof(int)); + xMemCpy(this->fglobaldoflist,pe->fglobaldoflist,this->fsize); } else{ @@ -265,5 +265,5 @@ /*}}}*/ /*FUNCTION ElementVector::SetValue{{{*/ -void ElementVector::SetValue(double scalar){ +void ElementVector::SetValue(IssmDouble scalar){ int i; Index: /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.h (revision 12472) @@ -22,5 +22,5 @@ int nrows; - double* values; + IssmDouble* values; //gset @@ -45,5 +45,5 @@ void CheckConsistency(void); void Init(ElementVector* pe); - void SetValue(double scalar); + void SetValue(IssmDouble scalar); /*}}}*/ }; Index: /issm/trunk-jpl/src/c/objects/Numerics/Matrix.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Numerics/Matrix.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Numerics/Matrix.cpp (revision 12472) @@ -45,10 +45,10 @@ #endif #ifdef _HAVE_ADOLC_ - this->amatrix=xNew(M*N); - #endif -} -/*}}}*/ -/*FUNCTION Matrix::Matrix(int M,int N,double sparsity){{{*/ -Matrix::Matrix(int M,int N,double sparsity){ + this->amatrix=xNew(M*N); + #endif +} +/*}}}*/ +/*FUNCTION Matrix::Matrix(int M,int N,IssmDouble sparsity){{{*/ +Matrix::Matrix(int M,int N,IssmDouble sparsity){ #ifdef _HAVE_PETSC_ @@ -58,10 +58,10 @@ #endif #ifdef _HAVE_ADOLC_ - this->amatrix=xNew; - #endif -} -/*}}}*/ -/*FUNCTION Matrix::Matrix(double* serial_mat, int M,int N,double sparsity){{{*/ -Matrix::Matrix(double* serial_mat, int M,int N,double sparsity){ + this->amatrix=xNew; + #endif +} +/*}}}*/ +/*FUNCTION Matrix::Matrix(IssmDouble* serial_mat, int M,int N,IssmDouble sparsity){{{*/ +Matrix::Matrix(IssmDouble* serial_mat, int M,int N,IssmDouble sparsity){ #ifdef _HAVE_PETSC_ @@ -85,5 +85,5 @@ #endif #ifdef _HAVE_ADOLC_ - this->amatrix=xNew(M*N); + this->amatrix=xNew(M*N); #endif } @@ -98,5 +98,5 @@ #endif #ifdef _HAVE_ADOLC_ - this->amatrix=xNew(M*N); + this->amatrix=xNew(M*N); #endif } @@ -111,5 +111,5 @@ #endif #ifdef _HAVE_ADOLC_ - xDelete(this->amatrix); + xDelete(this->amatrix); #endif } @@ -156,7 +156,7 @@ /*}}}*/ /*FUNCTION Matrix::Norm{{{*/ -double Matrix::Norm(NormMode norm_type){ - - double norm=0; +IssmDouble Matrix::Norm(NormMode norm_type){ + + IssmDouble norm=0; #ifdef _HAVE_PETSC_ _assert_(this->matrix); @@ -220,7 +220,7 @@ /*}}}*/ /*FUNCTION Matrix::ToSerial{{{*/ -double* Matrix::ToSerial(void){ - - double* output=NULL; +IssmDouble* Matrix::ToSerial(void){ + + IssmDouble* output=NULL; #ifdef _HAVE_PETSC_ @@ -233,5 +233,5 @@ /*}}}*/ /*FUNCTION Matrix::SetValues{{{*/ -void Matrix::SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode){ +void Matrix::SetValues(int m,int* idxm,int n,int* idxn,IssmDouble* values,InsMode mode){ #ifdef _HAVE_PETSC_ Index: /issm/trunk-jpl/src/c/objects/Numerics/Matrix.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Numerics/Matrix.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Numerics/Matrix.h (revision 12472) @@ -35,5 +35,5 @@ #endif #ifdef _HAVE_ADOLC_ - adouble* amatrix; + IssmDouble* amatrix; #endif @@ -41,6 +41,6 @@ Matrix(); Matrix(int M,int N); - Matrix(int M,int N,double sparsity); - Matrix(double* serial_mat,int M,int N,double sparsity); + Matrix(int M,int N,IssmDouble sparsity); + Matrix(IssmDouble* serial_mat,int M,int N,IssmDouble sparsity); Matrix(int M,int N,int connectivity,int numberofdofspernode); ~Matrix(); @@ -49,11 +49,11 @@ void Echo(void); void Assemble(void); - double Norm(NormMode norm_type); + IssmDouble Norm(NormMode norm_type); void GetSize(int* pM,int* pN); void GetLocalSize(int* pM,int* pN); void MatMult(Vector* X,Vector* AX); Matrix* Duplicate(void); - double* ToSerial(void); - void SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode); + IssmDouble* ToSerial(void); + void SetValues(int m,int* idxm,int n,int* idxn,IssmDouble* values,InsMode mode); void Convert(MatrixType type); /*}}}*/ Index: /issm/trunk-jpl/src/c/objects/Numerics/Vector.cpp =================================================================== --- /issm/trunk-jpl/src/c/objects/Numerics/Vector.cpp (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Numerics/Vector.cpp (revision 12472) @@ -44,10 +44,10 @@ #endif #ifdef _HAVE_ADOLC_ - this->avector=xNew(pM); - #endif -} -/*}}}*/ -/*FUNCTION Vector::Vector(double* serial_vec,int M){{{*/ -Vector::Vector(double* serial_vec,int M){ + this->avector=xNew(pM); + #endif +} +/*}}}*/ +/*FUNCTION Vector::Vector(IssmDouble* serial_vec,int M){{{*/ +Vector::Vector(IssmDouble* serial_vec,int M){ #ifdef _HAVE_PETSC_ @@ -65,5 +65,5 @@ #endif #ifdef _HAVE_ADOLC_ - this->avector=xNew(M); + this->avector=xNew(M); #endif } @@ -109,5 +109,5 @@ #endif #ifdef _HAVE_ADOLC_ - xDelete(this->avector); + xDelete(this->avector); #endif } @@ -146,5 +146,5 @@ /*}}}*/ /*FUNCTION Vector::SetValues{{{*/ -void Vector::SetValues(int ssize, int* list, double* values, InsMode mode){ +void Vector::SetValues(int ssize, int* list, IssmDouble* values, InsMode mode){ @@ -159,5 +159,5 @@ /*}}}*/ /*FUNCTION Vector::SetValue{{{*/ -void Vector::SetValue(int dof, double value, InsMode mode){ +void Vector::SetValue(int dof, IssmDouble value, InsMode mode){ #ifdef _HAVE_PETSC_ @@ -171,5 +171,5 @@ /*}}}*/ /*FUNCTION Vector::GetValue{{{*/ -void Vector::GetValue(double* pvalue,int dof){ +void Vector::GetValue(IssmDouble* pvalue,int dof){ #ifdef _HAVE_PETSC_ @@ -238,5 +238,5 @@ /*}}}*/ /*FUNCTION Vector::Set{{{*/ -void Vector::Set(double value){ +void Vector::Set(IssmDouble value){ #ifdef _HAVE_PETSC_ @@ -250,5 +250,5 @@ /*}}}*/ /*FUNCTION Vector::AXPY{{{*/ -void Vector::AXPY(Vector* X, double a){ +void Vector::AXPY(Vector* X, IssmDouble a){ #ifdef _HAVE_PETSC_ @@ -261,5 +261,5 @@ /*}}}*/ /*FUNCTION Vector::AYPX{{{*/ -void Vector::AYPX(Vector* X, double a){ +void Vector::AYPX(Vector* X, IssmDouble a){ #ifdef _HAVE_PETSC_ @@ -273,7 +273,7 @@ /*}}}*/ /*FUNCTION Vector::ToMPISerial{{{*/ -double* Vector::ToMPISerial(void){ - - double* vec_serial=NULL; +IssmDouble* Vector::ToMPISerial(void){ + + IssmDouble* vec_serial=NULL; #ifdef _HAVE_PETSC_ @@ -299,7 +299,7 @@ /*}}}*/ /*FUNCTION Vector::Norm{{{*/ -double Vector::Norm(NormMode norm_type){ - - double norm=0; +IssmDouble Vector::Norm(NormMode norm_type){ + + IssmDouble norm=0; #ifdef _HAVE_PETSC_ _assert_(this->vector); @@ -312,5 +312,5 @@ /*}}}*/ /*FUNCTION Vector::Scale{{{*/ -void Vector::Scale(double scale_factor){ +void Vector::Scale(IssmDouble scale_factor){ #ifdef _HAVE_PETSC_ @@ -323,7 +323,7 @@ /*}}}*/ /*FUNCTION Vector::Dot{{{*/ -double Vector::Dot(Vector* vector){ - - double dot; +IssmDouble Vector::Dot(Vector* vector){ + + IssmDouble dot; #ifdef _HAVE_PETSC_ _assert_(this->vector); Index: /issm/trunk-jpl/src/c/objects/Numerics/Vector.h =================================================================== --- /issm/trunk-jpl/src/c/objects/Numerics/Vector.h (revision 12471) +++ /issm/trunk-jpl/src/c/objects/Numerics/Vector.h (revision 12472) @@ -33,5 +33,5 @@ #endif #ifdef _HAVE_ADOLC_ - adouble* avector; + IssmDouble* avector; #endif @@ -39,5 +39,5 @@ Vector(); Vector(int M,bool fromlocalsize=false); - Vector(double* serial_vec,int pM); + Vector(IssmDouble* serial_vec,int pM); #ifdef _HAVE_PETSC_ Vector(Vec petsc_vec); @@ -50,21 +50,21 @@ /*Vector specific routines {{{*/ void Echo(void); - void AXPY(Vector *X, double a); - void AYPX(Vector *X, double a); + void AXPY(Vector *X, IssmDouble a); + void AYPX(Vector *X, IssmDouble a); void Assemble(void); void Copy(Vector *to); - double Dot(Vector *vector); + IssmDouble Dot(Vector *vector); Vector *Duplicate(void); - void GetValue(double *pvalue, int dof); + void GetValue(IssmDouble *pvalue, int dof); void GetSize(int *pM); void GetLocalSize(int *pM); bool IsEmpty(void); - double Norm(NormMode norm_type); + IssmDouble Norm(NormMode norm_type); void PointwiseDivide(Vector *x,Vector*y); - void Scale(double scale_factor); - void Set(double value); - void SetValues(int ssize, int *list, double*values, InsMode mode); - void SetValue(int dof, double value, InsMode mode); - double *ToMPISerial(void); + void Scale(IssmDouble scale_factor); + void Set(IssmDouble value); + void SetValues(int ssize, int *list, IssmDouble*values, InsMode mode); + void SetValue(int dof, IssmDouble value, InsMode mode); + IssmDouble *ToMPISerial(void); /*}}}*/ };