Changeset 12644

Timestamp:

07/17/12 14:48:07 (13 years ago)

Author:

seroussi

Message:

added transient spc for enthalpy (need to be improved to use iomodelcreateconstraints)

File:

• issm/trunk-jpl/src/c/modules/ModelProcessorx/Enthalpy/CreateConstraintsEnthalpy.cpp (modified) (3 diffs)

issm/trunk-jpl/src/c/modules/ModelProcessorx/Enthalpy/CreateConstraintsEnthalpy.cpp

@@ -9 +9 @@
 #include "../../../objects/objects.h"
 #include "../../../shared/shared.h"
+#include "../../../include/include.h"
 #include "../ModelProcessorx.h"
 
@@ -14 +15 @@
 
         /*Intermediary*/
-        int i;
-        int count;
+        int    i,j;
+        int    count;
         int    dim;
+        int    M,N;
         int    numberofvertices;
+        bool   spcpresent=false;
         double heatcapacity;
         double referencetemperature;
         
         /*Output*/
+        IssmDouble *spcvector  = NULL;
+        IssmDouble* times=NULL;
+        IssmDouble* values=NULL;
         Constraints* constraints = NULL;
-        SpcStatic*    spcstatic  = NULL;
 
         /*Fetch parameters: */
@@ -44 +49 @@
 
         /*Fetch data: */
-        double *spctemperature=NULL;
-        iomodel->FetchData(&spctemperature,NULL,NULL,ThermalSpctemperatureEnum);
+        iomodel->FetchData(&spcvector,&M,&N,ThermalSpctemperatureEnum);
 
-        /*Initialize counter*/
-        count=0;
+        //FIX ME: SHOULD USE IOMODELCREATECONSTRAINTS 
+        /*Transient or static?:*/
+        if(M==numberofvertices){
+                /*static: just create Constraints objects*/
+                count=0;
 
-        /*Create constraints from x,y,z: */
-        for (i=0;i<numberofvertices;i++){
-                /*keep only this partition's nodes:*/
-                if((iomodel->my_vertices[i])){
+                for (i=0;i<numberofvertices;i++){
+                        /*keep only this partition's nodes:*/
+                        if((iomodel->my_vertices[i])){
 
-                        if (!xIsNan<IssmDouble>(spctemperature[i])){
+                                if (!xIsNan<IssmDouble>(spcvector[i])){
 
-                                constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,heatcapacity*(spctemperature[i]-referencetemperature),EnthalpyAnalysisEnum));
-                                count++;
+                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,heatcapacity*(spcvector[i]-referencetemperature),EnthalpyAnalysisEnum));
+                                        count++;
 
+                                }
                         }
                 }
         }
+        else if (M==(numberofvertices+1)){
+                /*transient: create transient SpcTransient objects. Same logic, except we need to retrieve 
+                 * various times and values to initialize an SpcTransient object: */
+                count=0;
 
-        /*Free data: */
-        xDelete<double>(spctemperature);
-        
+                /*figure out times: */
+                times=xNew<IssmDouble>(N);
+                for(j=0;j<N;j++){
+                        times[j]=spcvector[(M-1)*N+j];
+                }
+                /*unit conversion: */
+                UnitConversion(times,N,ExtToIuEnum,TimeEnum);
+
+                /*Create constraints from x,y,z: */
+                for (i=0;i<numberofvertices;i++){
+                        
+                        /*keep only this partition's nodes:*/
+                        if((iomodel->my_vertices[i])){
+
+                                /*figure out times and values: */
+                                values=xNew<IssmDouble>(N);
+                                spcpresent=false;
+                                for(j=0;j<N;j++){
+                                        values[j]=heatcapacity*(spcvector[i*N+j]-referencetemperature);
+                                        if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
+                                }
+
+                                if(spcpresent){
+                                        constraints->AddObject(new SpcTransient(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,N,times,values,EnthalpyAnalysisEnum));
+                                        count++;
+                                }
+                                xDelete<IssmDouble>(values);
+                        }
+                }
+        }
+        else{
+                _error2_("Size of field " << EnumToStringx(ThermalSpctemperatureEnum) << " not supported");
+        }
+
+        /*Free ressources:*/
+        xDelete<IssmDouble>(spcvector);
+        xDelete<IssmDouble>(times);
+        xDelete<IssmDouble>(values);
+
         /*Assign output pointer: */
         *pconstraints=constraints;