Changeset 21908

Timestamp:

08/01/17 16:22:36 (8 years ago)

Author:

adhikari

Message:

CHG: extract components of moment of intertia tensor

Location:

issm/trunk-jpl/src/c

Files:

• analyses/EsaAnalysis.cpp (modified) (1 diff)
• analyses/SealevelriseAnalysis.cpp (modified) (1 diff)
• classes/FemModel.cpp (modified) (2 diffs)
• classes/FemModel.h (modified) (1 diff)
• cores/sealevelrise_core.cpp (modified) (1 diff)
• cores/sealevelrise_core_noneustatic.cpp (modified) (4 diffs)
• shared/Enum/EnumDefinitions.h (modified) (1 diff)
• shared/Enum/EnumToStringx.cpp (modified) (1 diff)
• shared/Enum/StringToEnumx.cpp (modified) (3 diffs)

issm/trunk-jpl/src/c/analyses/EsaAnalysis.cpp

@@ -80 +80 @@
                 U_elastic_local[i-lower_row]= (love_h[nl-1])/2.0/sin(alpha/2.0);
                 H_elastic_local[i-lower_row]= 0; 
-                IssmDouble Pn,Pn1,Pn2;
-                IssmDouble Pn_p,Pn_p1,Pn_p2;
+                //IssmDouble Pn,Pn1,Pn2;
+                //IssmDouble Pn_p,Pn_p1,Pn_p2;
+                IssmDouble Pn = 0.; 
+                IssmDouble Pn1 = 0.; 
+                IssmDouble Pn2 = 0.; 
+                IssmDouble Pn_p = 0.; 
+                IssmDouble Pn_p1 = 0.; 
+                IssmDouble Pn_p2 = 0.; 
+
                 for (int n=0;n<nl;n++) {
                         IssmDouble deltalove_U;

issm/trunk-jpl/src/c/analyses/SealevelriseAnalysis.cpp

@@ -127 +127 @@
                         U_elastic_local[i-lower_row]= (love_h[nl-1])/2.0/sin(alpha/2.0);
                         H_elastic_local[i-lower_row]= 0; 
-                        IssmDouble Pn,Pn1,Pn2;
-                        IssmDouble Pn_p,Pn_p1,Pn_p2;
+                        IssmDouble Pn = 0.; 
+                        IssmDouble Pn1 = 0.; 
+                        IssmDouble Pn2 = 0.; 
+                        IssmDouble Pn_p = 0.; 
+                        IssmDouble Pn_p1 = 0.; 
+                        IssmDouble Pn_p2 = 0.; 
+
                         for (int n=0;n<nl;n++) {
                                 IssmDouble deltalove_G;

issm/trunk-jpl/src/c/classes/FemModel.cpp

@@ -3920 +3920 @@
 }
 /*}}}*/
-void FemModel::SealevelriseRotationalFeedback(Vector<IssmDouble>* pSgo_rot, Vector<IssmDouble>* pSg_old, IssmDouble* latitude, IssmDouble* longitude, IssmDouble* radius){/*{{{*/
+void FemModel::SealevelriseRotationalFeedback(Vector<IssmDouble>* pSgo_rot, Vector<IssmDouble>* pSg_old, IssmDouble* pIxz, IssmDouble* pIyz, IssmDouble* pIzz, IssmDouble* latitude, IssmDouble* longitude, IssmDouble* radius){/*{{{*/
 
         /*serialized vectors:*/
@@ -3996 +3996 @@
         pSgo_rot->Assemble();
         
+        /*Assign output pointers:*/
+        *pIxz=moi_list[0];
+        *pIyz=moi_list[1];
+        *pIzz=moi_list[2];
+
         /*Free ressources:*/
         xDelete<IssmDouble>(Sg_old);

issm/trunk-jpl/src/c/classes/FemModel.h

@@ -137 +137 @@
                 void SealevelriseEustatic(Vector<IssmDouble>* pSgi, IssmDouble* peustatic, IssmDouble* latitude, IssmDouble* longitude, IssmDouble* radius);
                 void SealevelriseNonEustatic(Vector<IssmDouble>* pSgo, Vector<IssmDouble>* pSg_old, IssmDouble* latitude, IssmDouble* longitude, IssmDouble* radius,bool verboseconvolution);
-                void SealevelriseRotationalFeedback(Vector<IssmDouble>* pSgo_rot, Vector<IssmDouble>* pSg_old, IssmDouble* latitude, IssmDouble* longitude, IssmDouble* radius);
+                void SealevelriseRotationalFeedback(Vector<IssmDouble>* pSgo_rot, Vector<IssmDouble>* pSg_old, IssmDouble* pIxz, IssmDouble* pIyz, IssmDouble* pIzz, IssmDouble* latitude, IssmDouble* longitude, IssmDouble* radius);
                 void SealevelriseGeodetic(Vector<IssmDouble>* pUp, Vector<IssmDouble>* pNorth, Vector<IssmDouble>* pEast, Vector<IssmDouble>* pSg_old, IssmDouble* latitude, IssmDouble* longitude, IssmDouble* radius, IssmDouble* xx, IssmDouble* yy, IssmDouble* zz); 
                 IssmDouble SealevelriseOceanAverage(Vector<IssmDouble>* Sg);

issm/trunk-jpl/src/c/cores/sealevelrise_core.cpp

@@ -83 +83 @@
                 Sg=sealevelrise_core_noneustatic(femmodel,Sg_eustatic); //ocean loading tems  (2nd and 5th terms on the RHS of Farrel and Clark)
 
-                                
-        
                 /*compute other geodetic signatures, such as absolute sea level chagne, components of 3-D crustal motion: */
                 /*Initialize:*/

issm/trunk-jpl/src/c/cores/sealevelrise_core_noneustatic.cpp

@@ -37 +37 @@
         IssmDouble          *radius    = NULL;
         IssmDouble           eustatic;
+        IssmDouble                              Ixz; 
+        IssmDouble                              Iyz; 
+        IssmDouble                              Izz; 
 
         /*Recover some parameters: */
@@ -83 +86 @@
                         /*call rotational feedback  module: */
                         Sgo_rot = new Vector<IssmDouble>(gsize); Sgo_rot->Assemble();
-                        femmodel->SealevelriseRotationalFeedback(Sgo_rot,Sg_old,latitude,longitude,radius); 
+                        femmodel->SealevelriseRotationalFeedback(Sgo_rot,Sg_old,&Ixz,&Iyz,&Izz,latitude,longitude,radius); 
                         Sgo_rot->Assemble(); 
 
+                        /*save changes in inertia tensor as results: */
+                        femmodel->results->AddResult(new GenericExternalResult<IssmDouble>(femmodel->results->Size()+1,SealevelInertiaTensorXZEnum,Ixz));
+                        femmodel->results->AddResult(new GenericExternalResult<IssmDouble>(femmodel->results->Size()+1,SealevelInertiaTensorYZEnum,Iyz));
+                        femmodel->results->AddResult(new GenericExternalResult<IssmDouble>(femmodel->results->Size()+1,SealevelInertiaTensorZZEnum,Izz));
+                
                         Sgo->AXPY(Sgo_rot,1); 
+                        delete Sgo_rot; 
                 }
                 
@@ -93 +102 @@
         
                 /*Sg is the sum of the eustatic term, and the ocean terms: */
-                Sg_eustatic->Copy(Sg); Sg->AXPY(Sgo,1); 
+                Sg_eustatic->Copy(Sg); Sg->AXPY(Sgo,1); delete Sgo; 
                 Sg->Shift(-Sgo_oceanaverage);
 
@@ -112 +121 @@
                 /*some minor verbosing adjustment:*/
                 if(count>1)verboseconvolution=false;
-                
         }
         if(VerboseConvergence()) _printf0_("\n   total number of iterations: " << count-1 << "\n");

issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

@@ -795 +795 @@
         SealevelAbsoluteEnum,
         SealevelEustaticEnum,
+        SealevelInertiaTensorXZEnum,
+        SealevelInertiaTensorYZEnum,
+        SealevelInertiaTensorZZEnum,
         SealevelObsEnum,
         SealevelWeightsEnum,

issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

@@ -775 +775 @@
                 case SealevelAbsoluteEnum : return "SealevelAbsolute";
                 case SealevelEustaticEnum : return "SealevelEustatic";
+                case SealevelInertiaTensorXZEnum : return "SealevelInertiaTensorXZ";
+                case SealevelInertiaTensorYZEnum : return "SealevelInertiaTensorYZ";
+                case SealevelInertiaTensorZZEnum : return "SealevelInertiaTensorZZ";
                 case SealevelObsEnum : return "SealevelObs";
                 case SealevelWeightsEnum : return "SealevelWeights";

issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

@@ -793 +793 @@
               else if (strcmp(name,"SealevelAbsolute")==0) return SealevelAbsoluteEnum;
               else if (strcmp(name,"SealevelEustatic")==0) return SealevelEustaticEnum;
+              else if (strcmp(name,"SealevelInertiaTensorXZ")==0) return SealevelInertiaTensorXZEnum;
+              else if (strcmp(name,"SealevelInertiaTensorYZ")==0) return SealevelInertiaTensorYZEnum;
+              else if (strcmp(name,"SealevelInertiaTensorZZ")==0) return SealevelInertiaTensorZZEnum;
               else if (strcmp(name,"SealevelObs")==0) return SealevelObsEnum;
               else if (strcmp(name,"SealevelWeights")==0) return SealevelWeightsEnum;
@@ -872 +875 @@
               else if (strcmp(name,"Matenhancedice")==0) return MatenhancediceEnum;
               else if (strcmp(name,"Matestar")==0) return MatestarEnum;
-              else if (strcmp(name,"Matpar")==0) return MatparEnum;
-              else if (strcmp(name,"Node")==0) return NodeEnum;
-              else if (strcmp(name,"Numericalflux")==0) return NumericalfluxEnum;
          else stage=8;
    }
    if(stage==8){
-              if (strcmp(name,"Neumannflux")==0) return NeumannfluxEnum;
+              if (strcmp(name,"Matpar")==0) return MatparEnum;
+              else if (strcmp(name,"Node")==0) return NodeEnum;
+              else if (strcmp(name,"Numericalflux")==0) return NumericalfluxEnum;
+              else if (strcmp(name,"Neumannflux")==0) return NeumannfluxEnum;
               else if (strcmp(name,"Param")==0) return ParamEnum;
               else if (strcmp(name,"Moulin")==0) return MoulinEnum;
@@ -995 +998 @@
               else if (strcmp(name,"MaxAbsVz")==0) return MaxAbsVzEnum;
               else if (strcmp(name,"FloatingArea")==0) return FloatingAreaEnum;
-              else if (strcmp(name,"GroundedArea")==0) return GroundedAreaEnum;
-              else if (strcmp(name,"IceMass")==0) return IceMassEnum;
-              else if (strcmp(name,"IceVolume")==0) return IceVolumeEnum;
          else stage=9;
    }
    if(stage==9){
-              if (strcmp(name,"IceVolumeAboveFloatation")==0) return IceVolumeAboveFloatationEnum;
+              if (strcmp(name,"GroundedArea")==0) return GroundedAreaEnum;
+              else if (strcmp(name,"IceMass")==0) return IceMassEnum;
+              else if (strcmp(name,"IceVolume")==0) return IceVolumeEnum;
+              else if (strcmp(name,"IceVolumeAboveFloatation")==0) return IceVolumeAboveFloatationEnum;
               else if (strcmp(name,"TotalFloatingBmb")==0) return TotalFloatingBmbEnum;
               else if (strcmp(name,"TotalGroundedBmb")==0) return TotalGroundedBmbEnum;