Changeset 9110
- Timestamp:
- 07/25/11 13:12:44 (14 years ago)
- Location:
- issm/trunk/src
- Files:
-
- 9 added
- 1 deleted
- 16 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h
r9096 r9110 532 532 EpsvelEnum, 533 533 MeanvelEnum, 534 ParameteroutputEnum,534 Fake30Enum, 535 535 OutputfilenameEnum, 536 536 WaterfractionEnum, … … 542 542 QmuMassFluxNumProfilesEnum, 543 543 PartEnum, 544 MaxSteadystateIterationsEnum 544 MaxSteadystateIterationsEnum, 545 RequestedOutputsEnum, 546 BasalFrictionEnum, 547 ViscousHeatingEnum 545 548 }; 546 549 -
issm/trunk/src/c/Makefile.am
r9010 r9110 666 666 ./modules/Responsex/Responsex.h\ 667 667 ./modules/Responsex/Responsex.cpp\ 668 ./modules/RequestedOutputsx/RequestedOutputsx.h\ 669 ./modules/RequestedOutputsx/RequestedOutputsx.cpp\ 668 670 ./modules/Scotchx/Scotchx.cpp\ 669 671 ./modules/Scotchx/Scotchx.h\ … … 1341 1343 ./modules/Responsex/Responsex.h\ 1342 1344 ./modules/Responsex/Responsex.cpp\ 1345 ./modules/RequestedOutputsx/RequestedOutputsx.h\ 1346 ./modules/RequestedOutputsx/RequestedOutputsx.cpp\ 1343 1347 ./modules/Solverx/Solverx.cpp\ 1344 1348 ./modules/Solverx/Solverx.h\ -
issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp
r9096 r9110 473 473 case EpsvelEnum : return "Epsvel"; 474 474 case MeanvelEnum : return "Meanvel"; 475 case ParameteroutputEnum : return "Parameteroutput";475 case Fake30Enum : return "Fake30"; 476 476 case OutputfilenameEnum : return "Outputfilename"; 477 477 case WaterfractionEnum : return "Waterfraction"; … … 484 484 case PartEnum : return "Part"; 485 485 case MaxSteadystateIterationsEnum : return "MaxSteadystateIterations"; 486 case RequestedOutputsEnum : return "RequestedOutputs"; 487 case BasalFrictionEnum : return "BasalFriction"; 488 case ViscousHeatingEnum : return "ViscousHeating"; 486 489 default : return "unknown"; 487 490 -
issm/trunk/src/c/modules/ModelProcessorx/CreateParameters.cpp
r9013 r9110 18 18 Parameters* parameters = NULL; 19 19 20 char** parameteroutput=NULL;21 int num parameters=0;20 double* requestedoutputs=NULL; 21 int numoutputs=0; 22 22 23 23 if(*pparameters)return; //do not create parameters twice! … … 81 81 xfree((void**)&iomodel->riftinfo); 82 82 83 84 /*Requested output?*/ 85 IoModelFetchData(&requestedoutputs,&numoutputs,NULL,iomodel_handle,RequestedOutputsEnum); 86 parameters->AddObject(new IntVecParam(RequestedOutputsEnum,requestedoutputs,numoutputs)); 87 xfree((void**)&requestedoutputs); 83 88 84 /*Deal with parametr outputs: */85 IoModelFetchData(¶meteroutput,&numparameters,iomodel_handle,ParameteroutputEnum);86 if(numparameters) parameters->AddObject(new StringArrayParam(ParameteroutputEnum,parameteroutput,numparameters));87 88 89 /*Before returning, create parameters in case we are running Qmu or control types runs: */ 89 90 CreateParametersControl(¶meters,iomodel,iomodel_handle,solution_type,analysis_type); -
issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp
r9096 r9110 471 471 else if (strcmp(name,"Epsvel")==0) return EpsvelEnum; 472 472 else if (strcmp(name,"Meanvel")==0) return MeanvelEnum; 473 else if (strcmp(name," Parameteroutput")==0) return ParameteroutputEnum;473 else if (strcmp(name,"Fake30")==0) return Fake30Enum; 474 474 else if (strcmp(name,"Outputfilename")==0) return OutputfilenameEnum; 475 475 else if (strcmp(name,"Waterfraction")==0) return WaterfractionEnum; … … 482 482 else if (strcmp(name,"Part")==0) return PartEnum; 483 483 else if (strcmp(name,"MaxSteadystateIterations")==0) return MaxSteadystateIterationsEnum; 484 else if (strcmp(name,"RequestedOutputs")==0) return RequestedOutputsEnum; 485 else if (strcmp(name,"BasalFriction")==0) return BasalFrictionEnum; 486 else if (strcmp(name,"ViscousHeating")==0) return ViscousHeatingEnum; 484 487 else _error_("Enum %s not found",name); 485 488 -
issm/trunk/src/c/modules/modules.h
r9010 r9110 91 91 #include "./Reducevectorgtosx/Reducevectorgtosx.h" 92 92 #include "./Reducevectorgtofx/Reducevectorgtofx.h" 93 #include "./RequestedOutputsx/RequestedOutputsx.h" 93 94 #include "./Responsex/Responsex.h" 94 95 #include "./RheologyBbarx/RheologyBbarx.h" -
issm/trunk/src/c/objects/Elements/Element.h
r9012 r9110 67 67 virtual void ControlInputScaleGradient(int enum_type, double scale)=0; 68 68 virtual void ProcessResultsUnits(void)=0; 69 virtual void RequestedOutput(int output_enum,int step,double time)=0; 69 70 virtual void MinVel(double* pminvel, bool process_units)=0; 70 71 virtual void MaxVel(double* pmaxvel, bool process_units)=0; -
issm/trunk/src/c/objects/Elements/Penta.cpp
r9031 r9110 341 341 *(bed_normal+1)=-normal[1]/normal_norm; 342 342 *(bed_normal+2)=-normal[2]/normal_norm; 343 } 344 /*}}}*/ 345 /*FUNCTION Penta::BasalFrictionCreateInput {{{1*/ 346 void Penta::BasalFrictionCreateInput(void){ 347 348 /*Constants*/ 349 const int numdof=NUMVERTICES*NDOF1; 350 351 /*Intermediaries */ 352 int count,ig; 353 int drag_type; 354 double basalfriction[NUMVERTICES]={0,0,0,0,0,0}; 355 double alpha2,vx,vy; 356 Friction* friction=NULL; 357 GaussPenta* gauss=NULL; 358 359 360 /* Basal friction can only be found at the base of an ice sheet: */ 361 if (!IsOnBed() || IsOnShelf()){ 362 //empty friction: 363 this->inputs->AddInput(new PentaVertexInput(BasalFrictionEnum,&basalfriction[0])); 364 return; 365 } 366 367 /*Retrieve all inputs and parameters*/ 368 Input* vx_input=inputs->GetInput(VxEnum); _assert_(vx_input); 369 Input* vy_input=inputs->GetInput(VyEnum); _assert_(vy_input); 370 Input* vz_input=inputs->GetInput(VzEnum); _assert_(vz_input); 371 372 373 /*Build friction element, needed later: */ 374 inputs->GetParameterValue(&drag_type,DragTypeEnum); 375 if (drag_type!=2)_error_(" non-viscous friction not supported yet!"); 376 friction=new Friction("3d",inputs,matpar,DiagnosticHorizAnalysisEnum); 377 378 /* Start looping on the number of gauss 2d (nodes on the bedrock) */ 379 gauss=new GaussPenta(0,1,2,2); 380 count=0; 381 for(ig=gauss->begin();ig<gauss->end();ig++){ 382 383 gauss->GaussPoint(ig); 384 385 friction->GetAlpha2(&alpha2,gauss,VxEnum,VyEnum,VzEnum); 386 vx_input->GetParameterValue(&vx,gauss); 387 vy_input->GetParameterValue(&vy,gauss); 388 basalfriction[count]=alpha2*(pow(vx,2.0)+pow(vy,2.0)); 389 count++; 390 } 391 392 /*Create PentaVertex input, which will hold the basal friction:*/ 393 this->inputs->AddInput(new PentaVertexInput(BasalFrictionEnum,&basalfriction[0])); 394 395 /*Clean up and return*/ 396 delete gauss; 397 delete friction; 343 398 } 344 399 /*}}}*/ … … 5146 5201 if (enum_type==RheologyBbarEnum) input=this->matice->inputs->GetInput(RheologyBEnum); 5147 5202 else input=this->inputs->GetInput(enum_type); 5148 if (!input) _error_("Input %s not found in penta->inputs",EnumToStringx(enum_type)); 5203 //if (!input) _error_("Input %s not found in penta->inputs",EnumToStringx(enum_type)); why error out? if the requested input does not exist, we should still 5204 //try and output whatever we can instead of just failing. 5205 if(!input)return; 5149 5206 5150 5207 /*If we don't find it, no big deal, just don't do the transfer. Otherwise, build a new Result … … 7027 7084 /*Affect value to the reduced matrix */ 7028 7085 for(i=0;i<24;i++) *(Pe_reduced+i)=Pi[i]-Pright[i]; 7086 } 7087 /*}}}*/ 7088 /*FUNCTION Penta::RequestedOutput{{{1*/ 7089 void Penta::RequestedOutput(int output_enum,int step,double time){ 7090 7091 if(IsInput(output_enum)){ 7092 /*just transfer this input to results, and we are done: */ 7093 InputToResult(output_enum,step,time); 7094 } 7095 else{ 7096 /*this input does not exist, compute it, and then transfer to results: */ 7097 switch(output_enum){ 7098 case BasalFrictionEnum: 7099 7100 /*create input: */ 7101 BasalFrictionCreateInput(); 7102 7103 /*transfer to results :*/ 7104 InputToResult(output_enum,step,time); 7105 7106 /*erase input: */ 7107 inputs->DeleteInput(output_enum); 7108 break; 7109 case ViscousHeatingEnum: 7110 7111 /*create input: */ 7112 ViscousHeatingCreateInput(); 7113 7114 /*transfer to results :*/ 7115 InputToResult(output_enum,step,time); 7116 7117 /*erase input: */ 7118 inputs->DeleteInput(output_enum); 7119 break; 7120 7121 default: 7122 /*do nothing, no need to derail the computation because one of the outputs requested cannot be found: */ 7123 break; 7124 } 7125 } 7126 7029 7127 } 7030 7128 /*}}}*/ … … 7614 7712 } 7615 7713 /*}}}*/ 7714 /*FUNCTION Penta::ViscousHeatingCreateInput {{{1*/ 7715 void Penta::ViscousHeatingCreateInput(void){ 7716 7717 /*Constants*/ 7718 const int numdof=NUMVERTICES*NDOF1; 7719 7720 /*Intermediaries*/ 7721 int iv; 7722 double phi; 7723 double viscosity; 7724 double xyz_list[NUMVERTICES][3]; 7725 double epsilon[6]; 7726 double viscousheating[NUMVERTICES]={0,0,0,0,0,0}; 7727 GaussPenta *gauss=NULL; 7728 7729 /*Initialize Element vector*/ 7730 ElementVector* pe=new ElementVector(nodes,NUMVERTICES,this->parameters); 7731 7732 /*Retrieve all inputs and parameters*/ 7733 GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES); 7734 Input* vx_input=inputs->GetInput(VxEnum); _assert_(vx_input); 7735 Input* vy_input=inputs->GetInput(VyEnum); _assert_(vy_input); 7736 Input* vz_input=inputs->GetInput(VzEnum); _assert_(vz_input); 7737 7738 /*loop over vertices: */ 7739 gauss=new GaussPenta(); 7740 for (int iv=0;iv<NUMVERTICES;iv++){ 7741 gauss->GaussVertex(iv); 7742 7743 this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input); 7744 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]); 7745 GetPhi(&phi, &epsilon[0], viscosity); 7746 7747 viscousheating[iv]=phi; 7748 } 7749 7750 /*Create PentaVertex input, which will hold the basal friction:*/ 7751 this->inputs->AddInput(new PentaVertexInput(ViscousHeatingEnum,&viscousheating[0])); 7752 7753 /*Clean up and return*/ 7754 delete gauss; 7755 } 7756 /*}}}*/ -
issm/trunk/src/c/objects/Elements/Penta.h
r9012 r9110 75 75 /*Element virtual functions definitions: {{{1*/ 76 76 void AverageOntoPartition(Vec partition_contributions,Vec partition_areas,double* vertex_response,double* qmu_part); 77 void BasalFrictionCreateInput(void); 77 78 void ComputeBasalStress(Vec sigma_b); 78 79 void ComputeStrainRate(Vec eps); … … 120 121 void ShelfSync(); 121 122 double RheologyBbarAbsGradient(bool process_units,int weight_index); 123 void RequestedOutput(int output_enum,int step,double time); 122 124 double ThicknessAbsGradient(bool process_units,int weight_index); 123 125 void MigrateGroundingLine(); … … 143 145 int* GetHorizontalNeighboorSids(void); 144 146 double RheologyBbarx(void); 147 void ViscousHeatingCreateInput(void); 145 148 /*}}}*/ 146 149 /*Penta specific routines:{{{1*/ -
issm/trunk/src/c/objects/Elements/Tria.cpp
r9076 r9110 4679 4679 } 4680 4680 /*}}}*/ 4681 /*FUNCTION Tria::RequestedOutput{{{1*/ 4682 void Tria::RequestedOutput(int output_enum,int step,double time){ 4683 4684 if(IsInput(output_enum)){ 4685 /*just transfer this input to results, and we are done: */ 4686 InputToResult(output_enum,step,time); 4687 } 4688 else{ 4689 /*this input does not exist, compute it, and then transfer to results: */ 4690 switch(output_enum){ 4691 default: 4692 /*do nothing, no need to derail the computation because one of the outputs requested cannot be found: */ 4693 break; 4694 } 4695 } 4696 4697 } 4698 /*}}}*/ 4681 4699 /*FUNCTION Tria::SetClone {{{1*/ 4682 4700 void Tria::SetClone(int* minranks){ -
issm/trunk/src/c/objects/Elements/Tria.h
r9012 r9110 128 128 void MinVy(double* pminvy, bool process_units); 129 129 void MinVz(double* pminvz, bool process_units); 130 void RequestedOutput(int output_enum,int step,double time); 130 131 double RheologyBbarAbsGradient(bool process_units,int weight_index); 131 132 double ThicknessAbsMisfit( bool process_units,int weight_index); -
issm/trunk/src/c/objects/Params/IntVecParam.h
r8600 r9110 54 54 void GetParameterValue(int* pinteger){_error_("IntVec param of enum %i (%s) cannot return an integer",enum_type,EnumToStringx(enum_type));} 55 55 void GetParameterValue(int** pintarray,int* pM); 56 void GetParameterValue(int** pintarray,int* pM,int* pN){_error_("IntVec param of enum %i (%s) cannot return a string",enum_type,EnumToStringx(enum_type));}56 void GetParameterValue(int** pintarray,int* pM,int* pN){_error_("IntVec param of enum %i (%s) cannot return a matrix",enum_type,EnumToStringx(enum_type));} 57 57 void GetParameterValue(double* pdouble){_error_("IntVec param of enum %i (%s) cannot return a double",enum_type,EnumToStringx(enum_type));} 58 58 void GetParameterValue(char** pstring){_error_("IntVec param of enum %i (%s) cannot return a string",enum_type,EnumToStringx(enum_type));} -
issm/trunk/src/c/solutions/diagnostic_core.cpp
r9030 r9110 24 24 bool control_analysis; 25 25 int solution_type; 26 int numoutputs=0; 27 int* requested_outputs = NULL; 26 28 27 29 … … 34 36 femmodel->parameters->FindParam(&control_analysis,ControlAnalysisEnum); 35 37 femmodel->parameters->FindParam(&solution_type,SolutionTypeEnum); 38 femmodel->parameters->FindParam(&requested_outputs,&numoutputs,RequestedOutputsEnum); 36 39 37 40 /*for qmu analysis, reinitialize velocity so that fake sensitivities do not show up as a result of a different restart of the convergence at each trial.*/ … … 88 91 InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,PressureEnum); 89 92 if(dim==3) InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,VzEnum); 93 if(numoutputs)RequestedOutputsx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 90 94 } 95 96 /*Free ressources:*/ 97 xfree((void**)&requested_outputs); 98 91 99 } -
issm/trunk/src/c/solutions/steadystate_core.cpp
r9013 r9110 22 22 int max_steadystate_iterations; 23 23 bool control_analysis; 24 int numoutputs=0; 25 int* requested_outputs = NULL; 24 26 25 27 /* recover parameters:*/ … … 28 30 femmodel->parameters->FindParam(&solution_type,SolutionTypeEnum); 29 31 femmodel->parameters->FindParam(&max_steadystate_iterations,MaxSteadystateIterationsEnum); 32 femmodel->parameters->FindParam(&requested_outputs,&numoutputs,RequestedOutputsEnum); 30 33 31 34 /*intialize counters: */ … … 69 72 InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,TemperatureEnum); 70 73 InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,BasalMeltingRateEnum); 74 if(numoutputs)RequestedOutputsx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 71 75 } 76 77 /*Free ressources:*/ 78 xfree((void**)&requested_outputs); 79 72 80 } -
issm/trunk/src/m/classes/model.m
r9078 r9110 236 236 num_cm_responses = {0,true,'Integer'}; 237 237 %Output 238 parameteroutput = {{},true,'StringArray'};238 requested_outputs = {{},true,'DoubleMat',3}; 239 239 viscousheating = {NaN,false}; 240 240 pressure_elem = {NaN,false}; -
issm/trunk/src/m/model/ismodelselfconsistent.m
r9075 r9110 101 101 fields={'diagnostic_ref'}; 102 102 checksize(md,fields,[md.numberofnodes 6]); 103 if(size(md.requested_outputs,2)~=1), 104 message(['model ' md.name ' requested outputs should be a column vector']); 105 end 103 106 %}}} 104 107 %THICKNESS = SURFACE - BED {{{1 … … 259 262 end 260 263 if ~md.qmu_relax, 261 if md.eps_rel>1.1*10^-5,262 message(['model not consistent: for qmu analysis, eps_rel should be least than 10^-5, 10^-15 being a better value']);263 end264 %if md.eps_rel>1.1*10^-5, 265 % message(['model not consistent: for qmu analysis, eps_rel should be least than 10^-5, 10^-15 being a better value']); 266 %end 264 267 end 265 268 end
Note:
See TracChangeset
for help on using the changeset viewer.