Changeset 9880

Timestamp:

09/22/11 09:27:41 (13 years ago)

Author:

Mathieu Morlighem

Message:

Added Ice Volume as output

Location:

issm/trunk/src/c

Files:

• EnumDefinitions/EnumDefinitions.h (modified) (2 diffs)
• Makefile.am (modified) (2 diffs)
• modules/EnumToStringx/EnumToStringx.cpp (modified) (1 diff)
• modules/IceVolumex (added)
• modules/IceVolumex/IceVolumex.cpp (added)
• modules/IceVolumex/IceVolumex.h (added)
• modules/RequestedOutputsx/RequestedOutputsx.cpp (modified) (1 diff)
• modules/RequestedOutputsx/RequestedOutputsx.h (modified) (1 diff)
• modules/Responsex/Responsex.cpp (modified) (1 diff)
• modules/StringToEnumx/StringToEnumx.cpp (modified) (1 diff)
• modules/modules.h (modified) (1 diff)
• objects/Elements/Element.h (modified) (1 diff)
• objects/Elements/Penta.cpp (modified) (1 diff)
• objects/Elements/Penta.h (modified) (1 diff)
• objects/Elements/Tria.cpp (modified) (1 diff)
• objects/Elements/Tria.h (modified) (1 diff)
• solutions/diagnostic_core.cpp (modified) (1 diff)
• solutions/steadystate_core.cpp (modified) (1 diff)
• solutions/transient_core.cpp (modified) (1 diff)

issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h

@@ -475 +475 @@
         DragCoefficientAbsGradientEnum,
         TransientInputEnum,
-        /*Temporary*/
         OutputfilenameEnum,
         WaterfractionEnum,
@@ -485 +484 @@
         HydrologyWaterVyEnum,
         SpcDynamicEnum,
+        IceVolumeEnum,
         MaximumNumberOfEnums
 };

issm/trunk/src/c/Makefile.am

@@ -544 +544 @@
                                         ./modules/SurfaceAreax/SurfaceAreax.h\
                                         ./modules/SurfaceAreax/SurfaceAreax.cpp\
+                                        ./modules/IceVolumex/IceVolumex.h\
+                                        ./modules/IceVolumex/IceVolumex.cpp\
                                         ./modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.h\
                                         ./modules/SurfaceAbsVelMisfitx/SurfaceAbsVelMisfitx.cpp\
@@ -1241 +1243 @@
                                           ./modules/MinVzx/MinVzx.h\
                                           ./modules/MinVzx/MinVzx.cpp\
+                                          ./modules/IceVolumex/IceVolumex.h\
+                                          ./modules/IceVolumex/IceVolumex.cpp\
                                           ./modules/RheologyBbarx/RheologyBbarx.cpp\
                                           ./modules/RheologyBbarx/RheologyBbarx.h\

issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp

@@ -428 +428 @@
                 case HydrologyWaterVyEnum : return "HydrologyWaterVy";
                 case SpcDynamicEnum : return "SpcDynamic";
+                case IceVolumeEnum : return "IceVolume";
                 default : return "unknown";
 

issm/trunk/src/c/modules/RequestedOutputsx/RequestedOutputsx.cpp

@@ -9 +9 @@
 #include "../../EnumDefinitions/EnumDefinitions.h"
 
-void RequestedOutputsx(Elements* elements,Nodes* nodes,Vertices* vertices,Loads* loads,Materials* materials,Parameters* parameters, int* requested_outputs, int numoutputs, int step, double time){
+void RequestedOutputsx(Results* results,Elements* elements,Nodes* nodes,Vertices* vertices,Loads* loads,Materials* materials,Parameters* parameters, int* requested_outputs, int numoutputs, int step, double time){
 
-        int  i,o;
+        int  i,j;
         int  output_enum;
+        double output_value;
         Element* element=NULL;
 
         /*retrieve parameters: */
         if(numoutputs){
-                for(o=0;o<numoutputs;o++){
-                        output_enum=requested_outputs[o];
-                        /*create this output in the element inputs, and then transfer to results: */
-                        for(i=0;i<elements->Size();i++){
-                                element=(Element*)elements->GetObjectByOffset(i);
-                                element->RequestedOutput(output_enum,step,time);
+                for(i=0;i<numoutputs;i++){
+                        output_enum=requested_outputs[i];
+
+                        switch(output_enum){
+
+                                case IceVolumeEnum:
+                                        Responsex(&output_value,elements,nodes,vertices,loads,materials,parameters,"IceVolume",false,0);
+                                        results->AddObject(new DoubleExternalResult(results->Size()+1,IceVolumeEnum,output_value,step,time));
+                                        break;
+                                default:
+                                        /*create this output in the element inputs, and then transfer to results:*/
+                                        for(j=0;j<elements->Size();j++){
+                                                element=(Element*)elements->GetObjectByOffset(j);
+                                                element->RequestedOutput(output_enum,step,time);
+                                        }
+                                        break;
                         }
                 }

issm/trunk/src/c/modules/RequestedOutputsx/RequestedOutputsx.h

@@ -9 +9 @@
 
 /* local prototypes: */
-void RequestedOutputsx(Elements* elements,Nodes* nodes,Vertices* vertices,Loads* loads,Materials* materials,Parameters* parameters, int* requested_outputs, int numoutputs, int step=1, double time=0);
+void RequestedOutputsx(Results* results,Elements* elements,Nodes* nodes,Vertices* vertices,Loads* loads,Materials* materials,Parameters* parameters, int* requested_outputs, int numoutputs, int step=1, double time=0);
 
 #endif  /* _INPUTTORESULTX_H */

issm/trunk/src/c/modules/Responsex/Responsex.cpp

@@ -21 +21 @@
 
                 #ifdef _HAVE_RESPONSES_
+                case IceVolumeEnum:              IceVolumex(               responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
                 case MinVelEnum:                 MinVelx(                  responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;
                 case MaxVelEnum:                 MaxVelx(                  responses, elements,nodes, vertices, loads, materials, parameters,process_units); break;

issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp

@@ -426 +426 @@
         else if (strcmp(name,"HydrologyWaterVy")==0) return HydrologyWaterVyEnum;
         else if (strcmp(name,"SpcDynamic")==0) return SpcDynamicEnum;
+        else if (strcmp(name,"IceVolume")==0) return IceVolumeEnum;
         else _error_("Enum %s not found",name);
 

issm/trunk/src/c/modules/modules.h

@@ -33 +33 @@
 #include "./GroundingLineMigrationx/GroundingLineMigrationx.h"
 #include "./HoleFillerx/HoleFillerx.h"
+#include "./IceVolumex/IceVolumex.h"
 #include "./InputControlUpdatex/InputControlUpdatex.h"
 #include "./InputConvergencex/InputConvergencex.h"

issm/trunk/src/c/objects/Elements/Element.h

@@ -88 +88 @@
                 virtual double MassFlux(double* segment,bool process_units)=0;
                 virtual double RheologyBbarx(void)=0;
+                virtual double IceVolume(void)=0;
                 virtual int    NodalValue(double* pvalue, int index, int natureofdataenum,bool process_units)=0;
                 #endif

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

@@ -2667 +2667 @@
 
 #ifdef _HAVE_RESPONSES_
+/*FUNCTION Penta::IceVolume {{{1*/
+double Penta::IceVolume(void){
+
+        /*The volume of a troncated prism is base * 1/3 sum(length of edges)*/
+        double base,height;
+        double xyz_list[NUMVERTICES][3];
+
+        if(IsOnWater())return 0;
+
+        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
+
+        /*First calculate the area of the base (cross section triangle)
+         * http://en.wikipedia.org/wiki/Pentangle
+         * base = 1/2 abs((xA-xC)(yB-yA)-(xA-xB)(yC-yA))*/
+        base = 1./2.*fabs((xyz_list[0][0]-xyz_list[2][0])*(xyz_list[1][1]-xyz_list[0][1]) - (xyz_list[0][0]-xyz_list[1][0])*(xyz_list[2][1]-xyz_list[0][1]));
+
+        /*Now get the average height*/
+        height = 1./3.*((xyz_list[3][2]-xyz_list[0][2])+(xyz_list[4][2]-xyz_list[1][2])+(xyz_list[5][2]-xyz_list[2][2]));
+
+        /*Return: */
+        return base*height;
+}
+/*}}}*/
 /*FUNCTION Penta::RheologyBbarx{{{1*/
 double Penta::RheologyBbarx(void){

issm/trunk/src/c/objects/Elements/Penta.h

@@ -123 +123 @@
 
                  #ifdef _HAVE_RESPONSES_
+                double IceVolume(void);
                 void   MinVel(double* pminvel, bool process_units);
                 void   MinVx(double* pminvx, bool process_units);

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

@@ -2663 +2663 @@
 
 #ifdef _HAVE_RESPONSES_
+/*FUNCTION Tria::IceVolume {{{1*/
+double Tria::IceVolume(void){
+
+        /*The volume of a troncated prism is base * 1/3 sum(length of edges)*/
+        double base,surface,bed;
+        double xyz_list[NUMVERTICES][3];
+
+        if(IsOnWater())return 0;
+
+        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
+
+        /*First calculate the area of the base (cross section triangle)
+         * http://en.wikipedia.org/wiki/Triangle
+         * base = 1/2 abs((xA-xC)(yB-yA)-(xA-xB)(yC-yA))*/
+        base = 1./2. * fabs((xyz_list[0][0]-xyz_list[2][0])*(xyz_list[1][1]-xyz_list[0][1]) - (xyz_list[0][0]-xyz_list[1][0])*(xyz_list[2][1]-xyz_list[0][1]));
+
+        /*Now get the average height*/
+        Input* surface_input = inputs->GetInput(SurfaceEnum); _assert_(surface_input);
+        Input* bed_input     = inputs->GetInput(BedEnum);     _assert_(bed_input);
+        surface_input->GetParameterAverage(&surface);
+        bed_input->GetParameterAverage(&bed);
+
+        /*Return: */
+        return base*(surface-bed);
+}
+/*}}}*/
 /*FUNCTION Tria::MassFlux {{{1*/
 double Tria::MassFlux( double* segment,bool process_units){

issm/trunk/src/c/objects/Elements/Tria.h

@@ -124 +124 @@
 
                 #ifdef _HAVE_RESPONSES_
+                double IceVolume(void);
                 void   MinVel(double* pminvel, bool process_units);
                 void   MinVx(double* pminvx, bool process_units);

issm/trunk/src/c/solutions/diagnostic_core.cpp

@@ -92 +92 @@
                 InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,PressureEnum);
                 if(dim==3) InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,VzEnum);
-                RequestedOutputsx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,requested_outputs,numoutputs);
+                RequestedOutputsx(femmodel->results,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,requested_outputs,numoutputs);
         }
 

issm/trunk/src/c/solutions/steadystate_core.cpp

@@ -84 +84 @@
                 InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,TemperatureEnum);
                 InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,BasalforcingsMeltingRateEnum);
-                RequestedOutputsx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,requested_outputs,numoutputs);
+                RequestedOutputsx(femmodel->results,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,requested_outputs,numoutputs);
         }
 

issm/trunk/src/c/solutions/transient_core.cpp

@@ -122 +122 @@
                         InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,SurfaceforcingsMassBalanceEnum,step,time);
                         InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,MaskElementonfloatingiceEnum,step,time);
-                        RequestedOutputsx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,requested_outputs,numoutputs);
+                        RequestedOutputsx(femmodel->results,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,requested_outputs,numoutputs,step,time);
 
                         /*unload results*/