Changeset 548

Timestamp:

05/21/09 08:19:30 (16 years ago)

Author:

Mathieu Morlighem

Message:

missing parameters for thermal transient (initial temperature and pressure)

Location:

issm/trunk/src/c

Files:

• ModelProcessorx/Thermal/CreateParametersThermal.cpp (modified) (3 diffs)
• ProcessParamsx/ProcessParamsx.cpp (modified) (4 diffs)

issm/trunk/src/c/ModelProcessorx/Thermal/CreateParametersThermal.cpp

@@ -24 +24 @@
         double* vz=NULL;
         double* pressure=NULL;
+        double* temperature=NULL;
+        double* melting=NULL;
 
         /*recover parameters : */
@@ -30 +32 @@
         count=parameters->Size();
 
-        /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */
-        if (!model->ismacayealpattyn)return;
-
-        /*Get vx and vy: */
+        /*Get vx vy and vz: */
         ModelFetchData((void**)&vx,NULL,NULL,model_handle,"vx","Matrix","Mat");
         ModelFetchData((void**)&vy,NULL,NULL,model_handle,"vy","Matrix","Mat");
@@ -76 +75 @@
         xfree((void**)&pressure);
 
+        /* get initial temperature and melting if transient*/
+        if(strcmp(model->sub_analysis_type,"transient")==0){
+
+                /*Get melting and temperature: */
+                ModelFetchData((void**)&melting,NULL,NULL,model_handle,"melting","Matrix","Mat");
+                if(melting) for(i=0;i<model->numberofnodes;i++)melting[i]=melting[i]/model->yts;   //m/s instead of m/yr
+                ModelFetchData((void**)&temperature,NULL,NULL,model_handle,"temperature","Matrix","Mat");
+
+                count++;
+                param= new Param(count,"melting",DOUBLEVEC);
+                if(melting) param->SetDoubleVec(melting,model->numberofnodes);
+                else param->SetDoubleVec(melting,0);
+                parameters->AddObject(param);
+
+                count++;
+                param= new Param(count,"temperature",DOUBLEVEC);
+                if(temperature) param->SetDoubleVec(temperature,model->numberofnodes);
+                else throw ErrorException(__FUNCT__,exprintf("Missing initial temperature"));
+                parameters->AddObject(param);
+
+                /*Free melting and temperature: */
+                xfree((void**)&melting);
+                xfree((void**)&temperature);
+        }
         
         /*Assign output pointer: */
         *pparameters=parameters;
 }
-
-        

issm/trunk/src/c/ProcessParamsx/ProcessParamsx.cpp

@@ -24 +24 @@
         double* vy=NULL;
         double* vz=NULL;
+        
+        double* u_g=NULL;
+
+        /*thermal*/
         double* pressure=NULL;
-        
-        double* u_g=NULL;
+        double* temperature=NULL;
+        double* melting=NULL;
+
+        double* p_g=NULL;
+        double* t_g=NULL;
+        double* melting_g=NULL;
         
         /*control: */
@@ -38 +46 @@
 
         double* u_g_obs=NULL;
-        double* p_g=NULL;
-        
 
         /*diverse: */
         int     numberofnodes;
         int     analysis_type;
+        int     sub_analysis_type;
         int     count;
 
         parameters->FindParam((void*)&analysis_type,"analysis_type");
+        parameters->FindParam((void*)&sub_analysis_type,"sub_analysis_type");
         count=parameters->Size();
                 
@@ -162 +170 @@
                 param=(Param*)parameters->FindParamObject("pressure");
                 parameters->DeleteObject((Object*)param);
+
+                if(sub_analysis_type==TransientAnalysisEnum()){
+
+                        parameters->FindParam((void*)&temperature,"temperature");
+                        parameters->FindParam((void*)&melting,"melting");
+
+                        /*Now, from temperature and melting, build t_g and melting_g, correctly partitioned: */
+                        t_g=(double*)xcalloc(numberofnodes,sizeof(double));
+                        melting_g=(double*)xcalloc(numberofnodes,sizeof(double));
+
+                        for(i=0;i<numberofnodes;i++){
+                                t_g[(int)(partition[i])]=temperature[i];
+                                melting_g[(int)(partition[i])]=melting[i];
+                        }
+
+                        /*Now, create new parameter: */
+                        count++;
+                        param= new Param(count,"t_g",DOUBLEVEC);
+                        param->SetDoubleVec(t_g,numberofnodes);
+                        parameters->AddObject(param);
+
+                        count++;
+                        param= new Param(count,"melting_g",DOUBLEVEC);
+                        param->SetDoubleVec(melting_g,numberofnodes);
+                        parameters->AddObject(param);
+
+                        /*erase old parameter: */
+                        param=(Param*)parameters->FindParamObject("temperature");
+                        parameters->DeleteObject((Object*)param);
+
+                        param=(Param*)parameters->FindParamObject("melting");
+                        parameters->DeleteObject((Object*)param);
+                }
         }
 
@@ -174 +215 @@
         xfree((void**)&vy_obs);
         xfree((void**)&pressure);
+        xfree((void**)&temperature);
+        xfree((void**)&melting);
         xfree((void**)&control_parameter);
         xfree((void**)&u_g_obs);
         xfree((void**)&p_g); 
+        xfree((void**)&t_g); 
+        xfree((void**)&melting_g); 
 
 }