Changeset 9110


Ignore:
Timestamp:
07/25/11 13:12:44 (14 years ago)
Author:
Eric.Larour
Message:

Reworked parameteroutpus framework.
Now, we specify a requestedoutputs

Location:
issm/trunk/src
Files:
9 added
1 deleted
16 edited

Legend:

Unmodified
Added
Removed

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

    r9096 r9110  
    532532        EpsvelEnum,
    533533        MeanvelEnum,
    534         ParameteroutputEnum,
     534        Fake30Enum,
    535535        OutputfilenameEnum,
    536536        WaterfractionEnum,
    … …  
    542542        QmuMassFluxNumProfilesEnum,
    543543        PartEnum,
    544         MaxSteadystateIterationsEnum
     544        MaxSteadystateIterationsEnum,
     545        RequestedOutputsEnum,
     546        BasalFrictionEnum,
     547        ViscousHeatingEnum
    545548};
    546549

  • issm/trunk/src/c/Makefile.am

    r9010 r9110  
    666666                                        ./modules/Responsex/Responsex.h\
    667667                                        ./modules/Responsex/Responsex.cpp\
     668                                        ./modules/RequestedOutputsx/RequestedOutputsx.h\
     669                                        ./modules/RequestedOutputsx/RequestedOutputsx.cpp\
    668670                                        ./modules/Scotchx/Scotchx.cpp\
    669671                                        ./modules/Scotchx/Scotchx.h\
    … …  
    13411343                                        ./modules/Responsex/Responsex.h\
    13421344                                        ./modules/Responsex/Responsex.cpp\
     1345                                        ./modules/RequestedOutputsx/RequestedOutputsx.h\
     1346                                        ./modules/RequestedOutputsx/RequestedOutputsx.cpp\
    13431347                                        ./modules/Solverx/Solverx.cpp\
    13441348                                        ./modules/Solverx/Solverx.h\

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

    r9096 r9110  
    473473                case EpsvelEnum : return "Epsvel";
    474474                case MeanvelEnum : return "Meanvel";
    475                 case ParameteroutputEnum : return "Parameteroutput";
     475                case Fake30Enum : return "Fake30";
    476476                case OutputfilenameEnum : return "Outputfilename";
    477477                case WaterfractionEnum : return "Waterfraction";
    … …  
    484484                case PartEnum : return "Part";
    485485                case MaxSteadystateIterationsEnum : return "MaxSteadystateIterations";
     486                case RequestedOutputsEnum : return "RequestedOutputs";
     487                case BasalFrictionEnum : return "BasalFriction";
     488                case ViscousHeatingEnum : return "ViscousHeating";
    486489                default : return "unknown";
    487490

  • issm/trunk/src/c/modules/ModelProcessorx/CreateParameters.cpp

    r9013 r9110  
    1818        Parameters* parameters = NULL;
    1919
    20         char**   parameteroutput=NULL;
    21         int      numparameters=0;
     20        double*     requestedoutputs=NULL;
     21        int      numoutputs=0;
    2222       
    2323        if(*pparameters)return; //do not create parameters twice!
    … …  
    8181        xfree((void**)&iomodel->riftinfo);
    8282
     83
     84        /*Requested output?*/
     85        IoModelFetchData(&requestedoutputs,&numoutputs,NULL,iomodel_handle,RequestedOutputsEnum);
     86        parameters->AddObject(new IntVecParam(RequestedOutputsEnum,requestedoutputs,numoutputs));
     87        xfree((void**)&requestedoutputs);
    8388       
    84         /*Deal with parametr outputs: */
    85         IoModelFetchData(&parameteroutput,&numparameters,iomodel_handle,ParameteroutputEnum);
    86         if(numparameters) parameters->AddObject(new StringArrayParam(ParameteroutputEnum,parameteroutput,numparameters));
    87 
    8889        /*Before returning, create parameters in case we are running Qmu or control types runs: */
    8990        CreateParametersControl(&parameters,iomodel,iomodel_handle,solution_type,analysis_type);

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

    r9096 r9110  
    471471        else if (strcmp(name,"Epsvel")==0) return EpsvelEnum;
    472472        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;
    474474        else if (strcmp(name,"Outputfilename")==0) return OutputfilenameEnum;
    475475        else if (strcmp(name,"Waterfraction")==0) return WaterfractionEnum;
    … …  
    482482        else if (strcmp(name,"Part")==0) return PartEnum;
    483483        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;
    484487        else _error_("Enum %s not found",name);
    485488

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

    r9010 r9110  
    9191#include "./Reducevectorgtosx/Reducevectorgtosx.h"
    9292#include "./Reducevectorgtofx/Reducevectorgtofx.h"
     93#include "./RequestedOutputsx/RequestedOutputsx.h"
    9394#include "./Responsex/Responsex.h"
    9495#include "./RheologyBbarx/RheologyBbarx.h"

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

    r9012 r9110  
    6767                virtual void   ControlInputScaleGradient(int enum_type, double scale)=0;
    6868                virtual void   ProcessResultsUnits(void)=0;
     69                virtual void   RequestedOutput(int output_enum,int step,double time)=0;
    6970                virtual void   MinVel(double* pminvel, bool process_units)=0;
    7071                virtual void   MaxVel(double* pmaxvel, bool process_units)=0;

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

    r9031 r9110  
    341341        *(bed_normal+1)=-normal[1]/normal_norm;
    342342        *(bed_normal+2)=-normal[2]/normal_norm;
     343}
     344/*}}}*/
     345/*FUNCTION Penta::BasalFrictionCreateInput {{{1*/
     346void 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;
    343398}
    344399/*}}}*/
    … …  
    51465201        if (enum_type==RheologyBbarEnum) input=this->matice->inputs->GetInput(RheologyBEnum);
    51475202        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;
    51495206
    51505207        /*If we don't find it, no big deal, just don't do the transfer. Otherwise, build a new Result
    … …  
    70277084        /*Affect value to the reduced matrix */
    70287085        for(i=0;i<24;i++) *(Pe_reduced+i)=Pi[i]-Pright[i];
     7086}
     7087/*}}}*/
     7088/*FUNCTION Penta::RequestedOutput{{{1*/
     7089void 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
    70297127}
    70307128/*}}}*/
    … …  
    76147712}
    76157713/*}}}*/
     7714/*FUNCTION Penta::ViscousHeatingCreateInput {{{1*/
     7715void 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  
    7575                /*Element virtual functions definitions: {{{1*/
    7676                void   AverageOntoPartition(Vec partition_contributions,Vec partition_areas,double* vertex_response,double* qmu_part);
     77                void   BasalFrictionCreateInput(void);
    7778                void   ComputeBasalStress(Vec sigma_b);
    7879                void   ComputeStrainRate(Vec eps);
    … …  
    120121                void   ShelfSync();
    121122                double RheologyBbarAbsGradient(bool process_units,int weight_index);
     123                void   RequestedOutput(int output_enum,int step,double time);
    122124                double ThicknessAbsGradient(bool process_units,int weight_index);
    123125                void   MigrateGroundingLine();
    … …  
    143145                int*   GetHorizontalNeighboorSids(void);
    144146                double RheologyBbarx(void);
     147                void   ViscousHeatingCreateInput(void);
    145148                /*}}}*/
    146149                /*Penta specific routines:{{{1*/

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

    r9076 r9110  
    46794679}
    46804680/*}}}*/
     4681/*FUNCTION Tria::RequestedOutput{{{1*/
     4682void 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/*}}}*/
    46814699/*FUNCTION Tria::SetClone {{{1*/
    46824700void  Tria::SetClone(int* minranks){

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

    r9012 r9110  
    128128                void   MinVy(double* pminvy, bool process_units);
    129129                void   MinVz(double* pminvz, bool process_units);
     130                void   RequestedOutput(int output_enum,int step,double time);
    130131                double RheologyBbarAbsGradient(bool process_units,int weight_index);
    131132                double ThicknessAbsMisfit(     bool process_units,int weight_index);

  • issm/trunk/src/c/objects/Params/IntVecParam.h

    r8600 r9110  
    5454                void  GetParameterValue(int* pinteger){_error_("IntVec param of enum %i (%s) cannot return an integer",enum_type,EnumToStringx(enum_type));}
    5555                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));}
    5757                void  GetParameterValue(double* pdouble){_error_("IntVec param of enum %i (%s) cannot return a double",enum_type,EnumToStringx(enum_type));}
    5858                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  
    2424        bool   control_analysis;
    2525        int    solution_type;
     26        int    numoutputs=0;
     27        int*   requested_outputs    = NULL;
    2628
    2729
    … …  
    3436        femmodel->parameters->FindParam(&control_analysis,ControlAnalysisEnum);
    3537        femmodel->parameters->FindParam(&solution_type,SolutionTypeEnum);
     38        femmodel->parameters->FindParam(&requested_outputs,&numoutputs,RequestedOutputsEnum);
    3639
    3740        /*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.*/
    … …  
    8891                InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,PressureEnum);
    8992                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);
    9094        }
     95
     96        /*Free ressources:*/
     97        xfree((void**)&requested_outputs);
     98
    9199}

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

    r9013 r9110  
    2222        int max_steadystate_iterations;
    2323        bool control_analysis;
     24        int    numoutputs=0;
     25        int*   requested_outputs    = NULL;
    2426       
    2527        /* recover parameters:*/
    … …  
    2830        femmodel->parameters->FindParam(&solution_type,SolutionTypeEnum);
    2931        femmodel->parameters->FindParam(&max_steadystate_iterations,MaxSteadystateIterationsEnum);
     32        femmodel->parameters->FindParam(&requested_outputs,&numoutputs,RequestedOutputsEnum);
    3033
    3134        /*intialize counters: */
    … …  
    6972                InputToResultx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,TemperatureEnum);
    7073                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);
    7175        }
     76
     77        /*Free ressources:*/
     78        xfree((void**)&requested_outputs);
     79
    7280}

  • issm/trunk/src/m/classes/model.m

    r9078 r9110  
    236236                 num_cm_responses                    = {0,true,'Integer'};
    237237                 %Output
    238                  parameteroutput                     = {{},true,'StringArray'};
     238                 requested_outputs                    = {{},true,'DoubleMat',3};
    239239                 viscousheating                      = {NaN,false};
    240240                 pressure_elem                       = {NaN,false};

  • issm/trunk/src/m/model/ismodelselfconsistent.m

    r9075 r9110  
    101101fields={'diagnostic_ref'};
    102102checksize(md,fields,[md.numberofnodes 6]);
     103if(size(md.requested_outputs,2)~=1),
     104        message(['model ' md.name ' requested outputs should be a column vector']);
     105end
    103106%}}}
    104107%THICKNESS = SURFACE - BED {{{1
    … …  
    259262        end
    260263        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                 end
     264                %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
    264267        end
    265268end
Note: See TracChangeset for help on using the changeset viewer.