source: issm/trunk/src/c/ModelProcessorx/Control/CreateParametersControl.cpp@ 2268

Last change on this file since 2268 was 2268, checked in by Mathieu Morlighem, 15 years ago

Added regularization term in control methods. Misfit needs to be completed

File size: 6.1 KB
Line 
1/*!\file: CreateParametersControl.cpp
2 * \brief driver for creating parameters dataset, for control analysis.
3 */
4
5#undef __FUNCT__
6#define __FUNCT__ "CreateParameters"
7
8#include "../../DataSet/DataSet.h"
9#include "../../toolkits/toolkits.h"
10#include "../../EnumDefinitions/EnumDefinitions.h"
11#include "../../objects/objects.h"
12#include "../../shared/shared.h"
13#include "../IoModel.h"
14
15void CreateParametersControl(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle){
16
17 int i;
18
19 DataSet* parameters=NULL;
20 Param* param = NULL;
21 int count;
22 int analysis_type;
23 int numberofdofspernode;
24
25 double* fit=NULL;
26 double* cm_jump=NULL;
27 double* optscal=NULL;
28 double* maxiter=NULL;
29 double* control_parameter=NULL;
30 double* param_g=NULL;
31
32 double* vx_obs=NULL;
33 double* vy_obs=NULL;
34 double* u_g_obs=NULL;
35
36 double* vx=NULL;
37 double* vy=NULL;
38 double* vz=NULL;
39 double* u_g=NULL;
40
41 /*Get parameters: */
42 parameters=*pparameters;
43 count=parameters->Size();
44
45 //control analysis?
46 count++;
47 param= new Param(count,"control_analysis",INTEGER);
48 param->SetInteger(iomodel->control_analysis);
49 parameters->AddObject(param);
50
51 if(iomodel->control_analysis){
52 /*control_type: */
53 count++;
54 param= new Param(count,"control_type",STRING);
55 param->SetString(iomodel->control_type);
56 parameters->AddObject(param);
57
58 /*extrude_param: */
59 count++;
60 param= new Param(count,"extrude_param",DOUBLE);
61 if (strcmp(iomodel->control_type,"drag")==0) param->SetDouble(0);
62 else if (strcmp(iomodel->control_type,"B")==0) param->SetDouble(1);
63 else throw ErrorException(__FUNCT__,exprintf("control_type %s not supported yet!",iomodel->control_type));
64 parameters->AddObject(param);
65
66 /*control_steady: */
67 count++;
68 param= new Param(count,"control_steady",INTEGER);
69 param->SetInteger(0);
70 parameters->AddObject(param);
71
72 /*nsteps: */
73 count++;
74 param= new Param(count,"nsteps",INTEGER);
75 param->SetInteger(iomodel->nsteps);
76 parameters->AddObject(param);
77
78 /*tolx: */
79 count++;
80 param= new Param(count,"tolx",DOUBLE);
81 param->SetDouble(iomodel->tolx);
82 parameters->AddObject(param);
83
84 /*eps_cm: */
85 count++;
86 param= new Param(count,"eps_cm",DOUBLE);
87 param->SetDouble(iomodel->eps_cm);
88 parameters->AddObject(param);
89
90 /*cm_noisedampening: */
91 count++;
92 param= new Param(count,"cm_noisedampening",DOUBLE);
93 param->SetDouble(iomodel->cm_noisedampening);
94 parameters->AddObject(param);
95
96 /*mincontrolconstraint: */
97 count++;
98 param= new Param(count,"mincontrolconstraint",DOUBLE);
99 param->SetDouble(iomodel->mincontrolconstraint);
100 parameters->AddObject(param);
101
102 /*maxcontrolconstraint: */
103 count++;
104 param= new Param(count,"maxcontrolconstraint",DOUBLE);
105 param->SetDouble(iomodel->maxcontrolconstraint);
106 parameters->AddObject(param);
107
108 /*epsvel: */
109 count++;
110 param= new Param(count,"epsvel",DOUBLE);
111 param->SetDouble(iomodel->epsvel);
112 parameters->AddObject(param);
113
114 /*meanvel: */
115 count++;
116 param= new Param(count,"meanvel",DOUBLE);
117 param->SetDouble(iomodel->meanvel);
118 parameters->AddObject(param);
119
120 /*Now, recover fit, optscal and maxiter as vectors: */
121 IoModelFetchData((void**)&iomodel->fit,NULL,NULL,iomodel_handle,"fit","Matrix","Mat");
122 IoModelFetchData((void**)&iomodel->cm_jump,NULL,NULL,iomodel_handle,"cm_jump","Matrix","Mat");
123 IoModelFetchData((void**)&iomodel->optscal,NULL,NULL,iomodel_handle,"optscal","Matrix","Mat");
124 IoModelFetchData((void**)&iomodel->maxiter,NULL,NULL,iomodel_handle,"maxiter","Matrix","Mat");
125
126 count++;
127 param= new Param(count,"fit",DOUBLEVEC);
128 param->SetDoubleVec(iomodel->fit,iomodel->nsteps);
129 parameters->AddObject(param);
130
131 count++;
132 param= new Param(count,"cm_jump",DOUBLEVEC);
133 param->SetDoubleVec(iomodel->cm_jump,iomodel->nsteps);
134 parameters->AddObject(param);
135
136 count++;
137 param= new Param(count,"optscal",DOUBLEVEC);
138 param->SetDoubleVec(iomodel->optscal,iomodel->nsteps);
139 parameters->AddObject(param);
140
141 count++;
142 param= new Param(count,"maxiter",DOUBLEVEC);
143 param->SetDoubleVec(iomodel->maxiter,iomodel->nsteps);
144 parameters->AddObject(param);
145
146 xfree((void**)&iomodel->fit);
147 xfree((void**)&iomodel->cm_jump);
148 xfree((void**)&iomodel->optscal);
149 xfree((void**)&iomodel->maxiter);
150
151 /*Get vx, vx_obs, vy, vy_obs, and the parameter value: */
152 IoModelFetchData((void**)&vx,NULL,NULL,iomodel_handle,"vx","Matrix","Mat");
153 IoModelFetchData((void**)&vy,NULL,NULL,iomodel_handle,"vy","Matrix","Mat");
154 IoModelFetchData((void**)&vz,NULL,NULL,iomodel_handle,"vz","Matrix","Mat");
155 IoModelFetchData((void**)&vx_obs,NULL,NULL,iomodel_handle,"vx_obs","Matrix","Mat");
156 IoModelFetchData((void**)&vy_obs,NULL,NULL,iomodel_handle,"vy_obs","Matrix","Mat");
157 IoModelFetchData((void**)&control_parameter,NULL,NULL,iomodel_handle,iomodel->control_type,"Matrix","Mat");
158
159 u_g=(double*)xcalloc(iomodel->numberofnodes*3,sizeof(double));
160 if(vx)for(i=0;i<iomodel->numberofnodes;i++)u_g[3*i+0]=vx[i]/iomodel->yts;
161 if(vy)for(i=0;i<iomodel->numberofnodes;i++)u_g[3*i+1]=vy[i]/iomodel->yts;
162 if(vz)for(i=0;i<iomodel->numberofnodes;i++)u_g[3*i+2]=vz[i]/iomodel->yts;
163
164 count++;
165 param= new Param(count,"u_g",DOUBLEVEC);
166 param->SetDoubleVec(u_g,3*iomodel->numberofnodes,3);
167 parameters->AddObject(param);
168
169 u_g_obs=(double*)xcalloc(iomodel->numberofnodes*2,sizeof(double));
170 if(vx_obs)for(i=0;i<iomodel->numberofnodes;i++)u_g_obs[2*i+0]=vx_obs[i]/iomodel->yts;
171 if(vy_obs)for(i=0;i<iomodel->numberofnodes;i++)u_g_obs[2*i+1]=vy_obs[i]/iomodel->yts;
172
173 count++;
174 param= new Param(count,"u_g_obs",DOUBLEVEC);
175 param->SetDoubleVec(u_g_obs,2*iomodel->numberofnodes,2);
176 parameters->AddObject(param);
177
178 param_g=(double*)xcalloc(iomodel->numberofnodes,sizeof(double));
179 for(i=0;i<iomodel->numberofnodes;i++)param_g[i]=control_parameter[i];
180
181 count++;
182 param= new Param(count,"param_g",DOUBLEVEC);
183 param->SetDoubleVec(param_g,iomodel->numberofnodes,1);
184 parameters->AddObject(param);
185
186 xfree((void**)&vx);
187 xfree((void**)&vy);
188 xfree((void**)&vz);
189 xfree((void**)&u_g);
190 xfree((void**)&vx_obs);
191 xfree((void**)&vy_obs);
192 xfree((void**)&u_g_obs);
193 xfree((void**)&param_g);
194 xfree((void**)&control_parameter);
195 }
196
197 /*Assign output pointer: */
198 *pparameters=parameters;
199}
Note: See TracBrowser for help on using the repository browser.