Changeset 3612


Ignore:
Timestamp:
04/23/10 11:09:28 (15 years ago)
Author:
Eric.Larour
Message:

First new interface with new inputs

Location:
issm/trunk/src
Files:
43 added
9 deleted
67 edited

Legend:

Unmodified
Added
Removed

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

    r3529 r3612  
    1111/* local prototypes: */
    1212void    ComputeBasalStressx( Vec* pp_g,DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads, DataSet* materials,  DataSet* parameters,
    13                         ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                        int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _COMPUTEBASALSTRESSX_H */

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

    r3529 r3612  
    1111/* local prototypes: */
    1212void    ComputePressurex( Vec* pp_g,DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads, DataSet* materials,  DataSet* parameters,
    13                         ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                        int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _COMPUTEPRESSUREX_H */

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

    r3529 r3612  
    1111/* local prototypes: */
    1212void    ComputeStrainRatex(Vec* eps_g,DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads, DataSet* materials,  DataSet* parameters,
    13                         ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                        int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _COMPUTESTRAINRATEX_H */

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

    r3446 r3612  
    1111/* local prototypes: */
    1212void CostFunctionx( double* pJ, DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials, DataSet* parameters,
    13                         ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                        int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _MISFITX_H */

  • issm/trunk/src/c/DataSet/DataSet.cpp

    r3599 r3612  
    225225                        dataset->AddObject(tria);
    226226                }
     227                else if(enum_type==TriaVertexInputEnum){
     228                        TriaVertexInput* triavertexinput=NULL;
     229                        triavertexinput=new TriaVertexInput();
     230                        triavertexinput->Demarshall(&marshalled_dataset);
     231                        dataset->AddObject(triavertexinput);
     232                }
    227233                else if(enum_type==SingEnum){
    228234                        Sing* sing=NULL;
    … …  
    316322
    317323/*Specific methods*/
    318 /*FUNCTION DataSet::AddEinput{{{1*/
    319 int  DataSet::AddEinput(Einput* in_einput){
    320 
    321         /*First, go through dataset of inputs and check whether any input
    322          * with the same name is already in. If so, erase the corresponding
    323          * object before adding this new one: */
    324         vector<Object*>::iterator object;
    325         Einput* einput=NULL;
    326 
    327         for ( object=objects.begin() ; object < objects.end(); object++ ){
    328 
    329                 einput=(Einput*)(*object); //assume this is an inputs dataset
    330 
    331                 if (einput->EnumType()==in_einput->EnumType()){
    332                         this->DeleteObject(einput);
    333                         break;
    334                 }
    335         }
    336         this->AddObject(in_einput);
    337 }
    338 /*}}}*/
    339324/*FUNCTION DataSet::AddObject{{{1*/
    340325int  DataSet::AddObject(Object* object){

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

    r3599 r3612  
    1212#include <vector>
    1313#include "../toolkits/toolkits.h"
    14 #include "../objects/Einput.h"
    1514#include "../objects/Object.h"
    1615
    … …  
    4746                int   MarshallSize();
    4847                int   AddObject(Object* object);
    49                 int   AddEinput(Einput* in_einput);
    5048                int   DeleteObject(int id);
    5149                int   Size();

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

    r3446 r3612  
    1111/* local prototypes: */
    1212void Dux( Vec* pdu_g, DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials, DataSet* parameters,
    13                         ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                        int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _DUX_H */

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

    r3599 r3612  
    9999        InputEnum,
    100100        TriaVertexInputEnum,
     101        BoolInputEnum,
     102        IntInputEnum,
     103        DoubleInputEnum,
    101104        /*Params: */
    102105        ParamEnum,
    … …  
    124127        VxEnum,
    125128        VyEnum,
     129        VzEnum,
     130        VxAverageEnum,
     131        VyAverageEnum,
     132        VzAverageEnum,
     133        VxObsEnum,
     134        VyObsEnum,
     135        VzObsEnum,
     136        VxOldEnum,
     137        VyOldEnum,
     138        VzOldEnum,
     139        DhDtEnum,
     140        ThicknessEnum,
     141        SurfaceEnum,
     142        BedEnum,
     143        DragCoefficientEnum,
     144        DragPEnum,
     145        DragQEnum,
     146        DragTypeEnum,
     147        RheologyBEnum,
     148        RheologyNEnum,
     149        MeltingRateEnum,
     150        AccumulationRateEnum,
     151        GeothermalFluxEnum,
     152        ElementOnIceShelfEnum,
     153        ElementOnBedEnum,
     154        ElementOnWaterEnum,
     155        ElementOnSurfaceEnum,
     156        SurfaceAreaEnum,
     157        WeightsEnum,
     158        FitEnum,
     159        AdjointxEnum,
     160        AdjointyEnum,
     161        PressureEnum
    126162        /*}}}*/
    127163

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

    r3446 r3612  
    1111/* local prototypes: */
    1212void Gradjx( Vec* pgrad_g, int numberofnodes, DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials,  DataSet* parameters,
    13                         ParameterInputs* inputs,int analysis_type,int sub_analysis_type,char* control_type);
     13                        int analysis_type,int sub_analysis_type,char* control_type);
    1414
    1515#endif  /* _GRADJX_H */

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

    r3599 r3612  
    3030                                        ./objects/BamgOpts.cpp\
    3131                                        ./objects/Element.h\
    32                                         ./objects/ElementProperties.h\
    33                                         ./objects/ElementProperties.cpp\
    3432                                        ./objects/Model.h\
    3533                                        ./objects/Model.cpp\
    … …  
    6765                                        ./objects/TriaVertexInput.h\
    6866                                        ./objects/TriaVertexInput.cpp\
     67                                        ./objects/BoolInput.h\
     68                                        ./objects/BoolInput.cpp\
     69                                        ./objects/IntInput.h\
     70                                        ./objects/IntInput.cpp\
     71                                        ./objects/DoubleInput.h\
     72                                        ./objects/DoubleInput.cpp\
    6973                                        ./objects/Sing.h\
    7074                                        ./objects/Sing.cpp\
    … …  
    8286                                        ./objects/Input.cpp\
    8387                                        ./objects/Einput.h\
    84                                         ./objects/ParameterInputs.h\
    85                                         ./objects/ParameterInputs.cpp\
    8688                                        ./objects/Spc.cpp\
    8789                                        ./objects/Spc.h\
    … …  
    104106                                        ./DataSet/DataSet.cpp\
    105107                                        ./DataSet/DataSet.h\
     108                                        ./DataSet/Inputs.cpp\
     109                                        ./DataSet/Inputs.h\
    106110                                        ./shared/shared.h\
    107111                                        ./shared/Alloc/alloc.h\
    … …  
    200204                                        ./io/FetchParams.cpp\
    201205                                        ./io/FetchNodeSets.cpp\
    202                                         ./io/ParameterInputsInit.cpp\
    203206                                        ./io/pfopen.cpp\
    204207                                        ./io/pfclose.cpp\
    … …  
    272275                                        ./Gradjx/Gradjx.h\
    273276                                        ./Gradjx/Gradjx.cpp\
    274                                         ./UpdateFromInputsx/UpdateFromInputsx.h\
    275                                         ./UpdateFromInputsx/UpdateFromInputsx.cpp\
    276277                                        ./UpdateInputsx/UpdateInputsx.h\
    277278                                        ./UpdateInputsx/UpdateInputsx.cpp\
    … …  
    441442                                        ./objects/BamgOpts.cpp\
    442443                                        ./objects/Element.h\
    443                                         ./objects/ElementProperties.h\
    444                                         ./objects/ElementProperties.cpp\
    445444                                        ./objects/Model.h\
    446445                                        ./objects/Model.cpp\
    … …  
    478477                                        ./objects/TriaVertexInput.h\
    479478                                        ./objects/TriaVertexInput.cpp\
     479                                        ./objects/BoolInput.h\
     480                                        ./objects/BoolInput.cpp\
     481                                        ./objects/IntInput.h\
     482                                        ./objects/IntInput.cpp\
     483                                        ./objects/DoubleInput.h\
     484                                        ./objects/DoubleInput.cpp\
    480485                                        ./objects/Sing.h\
    481486                                        ./objects/Sing.cpp\
    … …  
    490495                                        ./objects/Numpar.h\
    491496                                        ./objects/Numpar.cpp\
    492                                         ./objects/Input.h\
    493                                         ./objects/Input.cpp\
    494497                                        ./objects/Einput.h\
    495                                         ./objects/ParameterInputs.h\
    496                                         ./objects/ParameterInputs.cpp\
    497498                                        ./objects/Spc.cpp\
    498499                                        ./objects/Spc.h\
    … …  
    515516                                        ./DataSet/DataSet.cpp\
    516517                                        ./DataSet/DataSet.h\
     518                                        ./DataSet/Inputs.cpp\
     519                                        ./DataSet/Inputs.h\
    517520                                        ./shared/shared.h\
    518521                                        ./shared/Threads/issm_threads.h\
    … …  
    608611                                        ./io/WriteParams.cpp\
    609612                                        ./io/FetchNodeSets.cpp\
    610                                         ./io/ParameterInputsInit.cpp\
    611613                                        ./io/pfopen.cpp\
    612614                                        ./io/pfclose.cpp\
    … …  
    680682                                        ./Gradjx/Gradjx.h\
    681683                                        ./Gradjx/Gradjx.cpp\
    682                                         ./UpdateFromInputsx/UpdateFromInputsx.h\
    683                                         ./UpdateFromInputsx/UpdateFromInputsx.cpp\
    684684                                        ./UpdateInputsx/UpdateInputsx.h\
    685685                                        ./UpdateInputsx/UpdateInputsx.cpp\

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

    r3446 r3612  
    1111/* local prototypes: */
    1212void Misfitx( double* pJ, DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials, DataSet* parameters,
    13                         ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                        int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _MISFITX_H */

  • issm/trunk/src/c/ModelProcessorx/DiagnosticHoriz/CreateElementsNodesAndMaterialsDiagnosticHoriz.cpp

    r3568 r3612  
    4343                IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface");
    4444                IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed");
    45                 IoModelFetchData(&iomodel->drag,NULL,NULL,iomodel_handle,"drag");
    46                 IoModelFetchData(&iomodel->p,NULL,NULL,iomodel_handle,"p");
    47                 IoModelFetchData(&iomodel->q,NULL,NULL,iomodel_handle,"q");
     45                IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient");
     46                IoModelFetchData(&iomodel->drag_p,NULL,NULL,iomodel_handle,"drag_p");
     47                IoModelFetchData(&iomodel->drag_q,NULL,NULL,iomodel_handle,"drag_q");
    4848                IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf");
    4949                IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater");
    5050                IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type");
    51                 IoModelFetchData(&iomodel->B,NULL,NULL,iomodel_handle,"B");
    52                 IoModelFetchData(&iomodel->n,NULL,NULL,iomodel_handle,"n");
    53                 IoModelFetchData(&iomodel->accumulation,NULL,NULL,iomodel_handle,"accumulation");
    54                 IoModelFetchData(&iomodel->melting,NULL,NULL,iomodel_handle,"melting");
     51                IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B");
     52                IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n");
    5553               
    5654                for (i=0;i<iomodel->numberofelements;i++){
    … …  
    7169                /*Free data : */
    7270                xfree((void**)&iomodel->elements);
     71                xfree((void**)&iomodel->elements_type);
    7372                xfree((void**)&iomodel->thickness);
    7473                xfree((void**)&iomodel->surface);
    7574                xfree((void**)&iomodel->bed);
    76                 xfree((void**)&iomodel->drag);
    77                 xfree((void**)&iomodel->p);
    78                 xfree((void**)&iomodel->q);
     75                xfree((void**)&iomodel->drag_coefficient);
     76                xfree((void**)&iomodel->drag_p);
     77                xfree((void**)&iomodel->drag_q);
     78                xfree((void**)&iomodel->rheology_B);
     79                xfree((void**)&iomodel->rheology_n);
    7980                xfree((void**)&iomodel->elementoniceshelf);
    8081                xfree((void**)&iomodel->elementonwater);
    81                 xfree((void**)&iomodel->B);
    82                 xfree((void**)&iomodel->n);
    83                 xfree((void**)&iomodel->elements_type);
    84 
     82               
    8583        }
    8684        else{ //        if (strcmp(meshtype,"2d")==0)
    … …  
    8886                /*Fetch data needed: */
    8987                IoModelFetchData(&iomodel->elements,NULL,NULL,iomodel_handle,"elements");
     88                IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type");
    9089                IoModelFetchData(&iomodel->thickness,NULL,NULL,iomodel_handle,"thickness");
    9190                IoModelFetchData(&iomodel->surface,NULL,NULL,iomodel_handle,"surface");
    9291                IoModelFetchData(&iomodel->bed,NULL,NULL,iomodel_handle,"bed");
    93                 IoModelFetchData(&iomodel->drag,NULL,NULL,iomodel_handle,"drag");
    94                 IoModelFetchData(&iomodel->p,NULL,NULL,iomodel_handle,"p");
    95                 IoModelFetchData(&iomodel->q,NULL,NULL,iomodel_handle,"q");
     92                IoModelFetchData(&iomodel->drag_coefficient,NULL,NULL,iomodel_handle,"drag_coefficient");
     93                IoModelFetchData(&iomodel->drag_p,NULL,NULL,iomodel_handle,"drag_p");
     94                IoModelFetchData(&iomodel->drag_q,NULL,NULL,iomodel_handle,"drag_q");
     95                IoModelFetchData(&iomodel->rheology_B,NULL,NULL,iomodel_handle,"rheology_B");
     96                IoModelFetchData(&iomodel->rheology_n,NULL,NULL,iomodel_handle,"rheology_n");
     97                IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater");
    9698                IoModelFetchData(&iomodel->elementoniceshelf,NULL,NULL,iomodel_handle,"elementoniceshelf");
    9799                IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed");
    98100                IoModelFetchData(&iomodel->elementonsurface,NULL,NULL,iomodel_handle,"elementonsurface");
    99                 IoModelFetchData(&iomodel->elements_type,NULL,NULL,iomodel_handle,"elements_type");
    100                 IoModelFetchData(&iomodel->B,NULL,NULL,iomodel_handle,"B");
    101                 IoModelFetchData(&iomodel->n,NULL,NULL,iomodel_handle,"n");
    102                 IoModelFetchData(&iomodel->accumulation,NULL,NULL,iomodel_handle,"accumulation");
    103                 IoModelFetchData(&iomodel->melting,NULL,NULL,iomodel_handle,"melting");
    104                 IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater");
    105                
     101       
    106102                for (i=0;i<iomodel->numberofelements;i++){
    107103                        if(iomodel->my_elements[i]){
    … …  
    120116                /*Free data: */
    121117                xfree((void**)&iomodel->elements);
     118                xfree((void**)&iomodel->elements_type);
    122119                xfree((void**)&iomodel->thickness);
    123120                xfree((void**)&iomodel->surface);
    124121                xfree((void**)&iomodel->bed);
    125                 xfree((void**)&iomodel->drag);
    126                 xfree((void**)&iomodel->p);
    127                 xfree((void**)&iomodel->q);
     122                xfree((void**)&iomodel->drag_coefficient);
     123                xfree((void**)&iomodel->drag_p);
     124                xfree((void**)&iomodel->drag_q);
     125                xfree((void**)&iomodel->rheology_n);
     126                xfree((void**)&iomodel->rheology_B);
    128127                xfree((void**)&iomodel->elementoniceshelf);
    129128                xfree((void**)&iomodel->elementonbed);
    130129                xfree((void**)&iomodel->elementonsurface);
    131                 xfree((void**)&iomodel->elements_type);
    132                 xfree((void**)&iomodel->n);
    133                 xfree((void**)&iomodel->B);
    134                 xfree((void**)&iomodel->accumulation);
    135                 xfree((void**)&iomodel->melting);
    136130                xfree((void**)&iomodel->elementonwater);
     131
    137132
    138133        } //if (strcmp(meshtype,"2d")==0)

  • issm/trunk/src/c/ModelProcessorx/IoModel.cpp

    r3582 r3612  
    8585
    8686        iomodel->drag_type=0;
    87         iomodel->drag=NULL;
    88         iomodel->p=NULL;
    89         iomodel->q=NULL;
     87        iomodel->drag_cofficient=NULL;
     88        iomodel->drag_p=NULL;
     89        iomodel->drag_q=NULL;
    9090       
    9191       
    … …  
    102102        iomodel->rho_ice=0;
    103103        iomodel->g=0;
    104         iomodel->n=NULL;
    105         iomodel->B=NULL;
     104        iomodel->rheology_n=NULL;
     105        iomodel->rheology_B=NULL;
    106106
    107107        /*!control methods: */
    … …  
    170170
    171171        /*!basal: */
    172         iomodel->accumulation=NULL;
     172        iomodel->accumulation_rate=NULL;
    173173        iomodel->dhdt=NULL;
    174174       
    … …  
    240240        xfree((void**)&iomodel->pressure);
    241241        xfree((void**)&iomodel->temperature);
    242         xfree((void**)&iomodel->drag);
    243         xfree((void**)&iomodel->p);
    244         xfree((void**)&iomodel->q);
     242        xfree((void**)&iomodel->drag_coefficient);
     243        xfree((void**)&iomodel->drag_p);
     244        xfree((void**)&iomodel->drag_q);
    245245        xfree((void**)&iomodel->elementoniceshelf);
    246246        xfree((void**)&iomodel->elementonwater);
    … …  
    253253        xfree((void**)&iomodel->edges);
    254254        xfree((void**)&iomodel->geothermalflux);
    255         xfree((void**)&iomodel->melting);
    256         xfree((void**)&iomodel->accumulation);
     255        xfree((void**)&iomodel->melting_rate);
     256        xfree((void**)&iomodel->accumulation_rate);
    257257        xfree((void**)&iomodel->dhdt);
    258         xfree((void**)&iomodel->B);
    259         xfree((void**)&iomodel->n);
     258        xfree((void**)&iomodel->rheology_B);
     259        xfree((void**)&iomodel->rheology_n);
    260260        xfree((void**)&iomodel->fit);
    261261        xfree((void**)&iomodel->weights);

  • issm/trunk/src/c/ModelProcessorx/IoModel.h

    r3588 r3612  
    8484        /*friction: */
    8585        int     drag_type;
    86         double* drag;
    87         double* p;
    88         double* q;
     86        double* drag_coefficient;
     87        double* drag_p;
     88        double* drag_q;
    8989
    9090        /*boundary conditions: */
    … …  
    101101        double  rho_water,rho_ice;
    102102        double  g;
    103         double* B;
    104         double* n;
     103        double* rheology_B;
     104        double* rheology_n;
    105105
    106106        /*numerical parameters: */
    … …  
    170170       
    171171        /*basal: */
    172         double*  melting;
    173         double*  accumulation;
     172        double*  melting_rate;
     173        double*  accumulation_rate;
    174174        double*  dhdt;
    175175

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

    r3446 r3612  
    1111/* local prototypes: */
    1212void PenaltyConstraintsx(int* pconverged, int* pnum_unstable_constraints, DataSet* elements,DataSet* nodes, DataSet* vertices,
    13                 DataSet* loads,DataSet* materials,  DataSet* parameters,ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                DataSet* loads,DataSet* materials,  DataSet* parameters,int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _PENALTYCONSTRAINTSX_H */

  • issm/trunk/src/c/PenaltySystemMatricesx/PenaltySystemMatricesx.cpp

    r3595 r3612  
    1111
    1212void PenaltySystemMatricesx(Mat Kgg, Vec pg,double* pkmax,DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads,DataSet* materials, DataSet* parameters,
    13                 int kflag,int pflag,ParameterInputs* inputs,int analysis_type,int sub_analysis_type){
     13                int kflag,int pflag,int analysis_type,int sub_analysis_type){
    1414       
    1515        int i;

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

    r3446 r3612  
    1111/* local prototypes: */
    1212void PenaltySystemMatricesx(Mat Kgg, Vec pg,double* pkmax, DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads,DataSet* materials, DataSet* parameters,
    13                 int kflag,int pflag,ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                int kflag,int pflag,int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _PENALTYSYSTEMMATRICESX_H */

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

    r3446 r3612  
    1515void SpawnCoreSerial(double* responses, int numresponses, double* variables, char** variables_descriptors,int numvariables, mxArray* model,mxArray* inputs,int analysis_type,int sub_analysis_type,int counter);
    1616#else
    17 void Qmux(Model* model,ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
    18 void SpawnCoreParallel(double* responses, int numresponses, double* variables, char** variables_descriptors,int numvariables, Model* model,ParameterInputs* inputs,int analysis_type,int sub_analysis_type,int counter);
    19 void DakotaResponses(double* responses,char** responses_descriptors,int numresponses,Model* model, DataSet* results,DataSet* processed_results,ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     17void Qmux(Model* model,int analysis_type,int sub_analysis_type);
     18void SpawnCoreParallel(double* responses, int numresponses, double* variables, char** variables_descriptors,int numvariables, Model* model,int analysis_type,int sub_analysis_type,int counter);
     19void DakotaResponses(double* responses,char** responses_descriptors,int numresponses,Model* model, DataSet* results,DataSet* processed_results,int analysis_type,int sub_analysis_type);
    2020#endif
    2121

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

    r3446 r3612  
    1111/* local prototypes: */
    1212void SurfaceAreax( double* pS, DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials, DataSet* parameters,
    13                         ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                        int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _SURFACEAREAX_H */

  • issm/trunk/src/c/SystemMatricesx/SystemMatricesx.cpp

    r3483 r3612  
    1111
    1212void SystemMatricesx(Mat* pKgg, Vec* ppg,DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads,DataSet* materials, DataSet* parameters,
    13                 int kflag,int pflag,int connectivity,int numberofdofspernode,ParameterInputs* inputs,int analysis_type,int sub_analysis_type){
     13                int kflag,int pflag,int connectivity,int numberofdofspernode,int analysis_type,int sub_analysis_type){
    1414       
    1515        int i;
    … …  
    3030        loads->Configure(elements, loads, nodes,vertices, materials,parameters);
    3131        parameters->Configure(elements,loads, nodes,vertices, materials,parameters);
     32
    3233
    3334        /*Get size of matrix: */

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

    r3446 r3612  
    1111/* local prototypes: */
    1212void SystemMatricesx(Mat* pKgg, Vec* ppg,DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads,DataSet* materials, DataSet* parameters,
    13                 int kflag,int pflag,int connectivity,int numberofdofspernode,ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     13                int kflag,int pflag,int connectivity,int numberofdofspernode,int analysis_type,int sub_analysis_type);
    1414
    1515#endif  /* _SYSTEMMATRICESX_H */

  • issm/trunk/src/c/UpdateInputsx/UpdateInputsx.cpp

    r3599 r3612  
    3131        elements->UpdateInputs(serial_solution,analysis_type,sub_analysis_type);
    3232
    33         elements->Echo();
    34 
    3533        /*Free ressources:*/
    3634        xfree((void**)&serial_solution);

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

    r3446 r3612  
    77
    88#include "../DataSet/DataSet.h"
    9 #include "../objects/ParameterInputs.h"
    109
    1110/* local prototypes: */

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

    r3599 r3612  
    3939#include "./ConfigureObjectsx/ConfigureObjectsx.h"
    4040#include "./SystemMatricesx/SystemMatricesx.h"
    41 #include "./UpdateFromInputsx/UpdateFromInputsx.h"
    4241#include "./UpdateInputsx/UpdateInputsx.h"
    4342#include "./UpdateGeometryx/UpdateGeometryx.h"

  • issm/trunk/src/c/objects/Beam.cpp

    r3599 r3612  
    243243}
    244244/*}}}*/
    245 /*FUNCTION Beam::UpdateFromInputs{{{1*/
    246 void  Beam::UpdateFromInputs(void* vinputs){
    247 
    248         int     dofs[1]={0};
    249         double  temperature_list[2];
    250         double  temperature_average;
    251         double  B_list[2];
    252         double  B_average;
    253 
    254         /*dynamic objects pointed to by hooks: */
    255         Node**  nodes=NULL;
    256         Matpar* matpar=NULL;
    257         Matice* matice=NULL;
    258         Numpar* numpar=NULL;
    259 
    260         ParameterInputs* inputs=NULL;
    261 
    262         /*recover pointers: */
    263         inputs=(ParameterInputs*)vinputs;
    264 
    265         /*recover objects from hooks: */
    266         nodes=(Node**)hnodes.deliverp();
    267         matpar=(Matpar*)hmatpar.delivers();
    268         matice=(Matice*)hmatice.delivers();
    269         numpar=(Numpar*)hnumpar.delivers();
    270 
    271         /*Update internal data if inputs holds new values: */
    272         //if (id==1) printf("WARNING if QMU: no hydrostatic equilibrium is applied here (conflict with prognostic, which does not have matpar)\n");
    273         //For now
    274         if(this->properties.h) inputs->Recover("thickness",&this->properties.h[0],1,dofs,2,(void**)nodes);
    275         //Later
    276         /*
    277                 if(inputs->Recover("thickness",&new_h[0],1,dofs,2,(void**)nodes)){
    278         //density, needed later:
    279         double di=(this->matpar->GetRhoIce()/this->matpar->GetRhoWater());
    280         //Go through grids:
    281         for (i=0;i<2;i++){
    282         if(nodes[i]->IsOnShelf()){
    283         this->b[i]=this->b[i]-di*(new_h[i]-h[i]); //hydrostatic equilibrium;
    284         }
    285         this->s[i]=this->b[i]+new_h[i];
    286         this->h[i]=new_h[i];
    287         }
    288         }
    289         */
    290         if (this->properties.s)              inputs->Recover("surface",&this->properties.s[0],1,dofs,2,(void**)nodes);
    291         if (this->properties.b)              inputs->Recover("bed",&this->properties.b[0],1,dofs,2,(void**)nodes);
    292         if (this->properties.k)              inputs->Recover("drag",&this->properties.k[0],1,dofs,2,(void**)nodes);
    293         if (this->properties.melting)        inputs->Recover("melting",&this->properties.melting[0],1,dofs,2,(void**)nodes);
    294         if (this->properties.accumulation)   inputs->Recover("accumulation",&this->properties.accumulation[0],1,dofs,2,(void**)nodes);
    295         if (this->properties.geothermalflux) inputs->Recover("geothermalflux",&this->properties.geothermalflux[0],1,dofs,2,(void**)nodes);
    296        
    297         //Update material if necessary
    298         if(inputs->Recover("temperature_average",&temperature_list[0],1,dofs,2,(void**)nodes)){
    299                 temperature_average=(temperature_list[0]+temperature_list[1])/2.0;
    300                 B_average=Paterson(temperature_average);
    301                 matice->SetB(B_average);
    302         }
    303         if(inputs->Recover("temperature",&temperature_list[0],1,dofs,2,(void**)nodes)){
    304                 temperature_average=(temperature_list[0]+temperature_list[1])/2.0;
    305                 B_average=Paterson(temperature_average);
    306                 matice->SetB(B_average);
    307         }
    308         if(inputs->Recover("B",&B_list[0],1,dofs,2,(void**)nodes)){
    309                 B_average=(B_list[0]+B_list[1])/2.0;
    310                 matice->SetB(B_average);
    311         }
    312 
    313 }
    314 /*}}}*/
    315245
    316246/*Object functions*/
    317247/*FUNCTION Beam::ComputeBasalStress{{{1*/
    318 void  Beam::ComputeBasalStress(Vec eps,void* inputs,int analysis_type,int sub_analysis_type){
     248void  Beam::ComputeBasalStress(Vec eps,int analysis_type,int sub_analysis_type){
    319249
    320250        ISSMERROR("Not implemented yet");
    … …  
    323253/*}}}*/
    324254/*FUNCTION Beam::ComputePressure{{{1*/
    325 void  Beam::ComputePressure(Vec p_g,void* inputs,int analysis_type,int sub_analysis_type){
     255void  Beam::ComputePressure(Vec p_g,int analysis_type,int sub_analysis_type){
    326256
    327257        int i;
    … …  
    366296/*}}}*/
    367297/*FUNCTION Beam::ComputeStrainRate{{{1*/
    368 void  Beam::ComputeStrainRate(Vec eps,void* inputs,int analysis_type,int sub_analysis_type){
     298void  Beam::ComputeStrainRate(Vec eps,int analysis_type,int sub_analysis_type){
    369299
    370300        ISSMERROR("Not implemented yet");
    … …  
    379309/*FUNCTION Beam::CreateKMatrix{{{1*/
    380310
    381 void  Beam::CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     311void  Beam::CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type){
    382312
    383313        /*Just branch to the correct element stiffness matrix generator, according to the type of analysis we are carrying out: */
    … …  
    399329/*FUNCTION Beam::CreateKMatrixDiagnosticHutter{{{1*/
    400330
    401 void  Beam::CreateKMatrixDiagnosticHutter(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     331void  Beam::CreateKMatrixDiagnosticHutter(Mat Kgg,int analysis_type,int sub_analysis_type){
    402332       
    403333       
    … …  
    432362/*}}}*/
    433363/*FUNCTION Beam::CreatePVector{{{1*/
    434 void  Beam::CreatePVector(Vec pg,void* inputs,int analysis_type,int sub_analysis_type){
     364void  Beam::CreatePVector(Vec pg,int analysis_type,int sub_analysis_type){
    435365       
    436366        /*Just branch to the correct load generator, according to the type of analysis we are carrying out: */
    … …  
    450380/*FUNCTION Beam::CreatePVectorDiagnosticHutter{{{1*/
    451381
    452 void Beam::CreatePVectorDiagnosticHutter( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type){
     382void Beam::CreatePVectorDiagnosticHutter( Vec pg, int analysis_type,int sub_analysis_type){
    453383
    454384        int i,j,k;
    … …  
    699629/*}}}*/
    700630/*FUNCTION Beam::GetOnBed{{{1*/
    701 int   Beam::GetOnBed(){
     631bool   Beam::GetOnBed(){
    702632        ISSMERROR(" not supported yet!");
    703633}
    … …  
    714644/*}}}*/
    715645/*FUNCTION Beam::GetShelf{{{1*/
    716 int   Beam::GetShelf(){
     646bool   Beam::GetShelf(){
    717647        ISSMERROR(" not supported yet!");
    718648}

  • issm/trunk/src/c/objects/Beam.h

    r3599 r3612  
    1111#include "./Matice.h"
    1212#include "./Matpar.h"
    13 #include "./ParameterInputs.h"
    14 #include "./ElementProperties.h"
    1513#include "../ModelProcessorx/IoModel.h"
    1614#include "./Hook.h"
    … …  
    2018class Element;
    2119class Hook;
    22 class ElementProperties;
    2320
    2421class Beam: public Element{
    … …  
    3431                Hook hnumpar; //hook to 1 numpar
    3532
    36                 ElementProperties properties;
     33                Inputs* inputs;
    3734       
    3835        public:
    … …  
    4037                /*constructors, destructors: {{{1*/
    4138                Beam();
    42                 Beam(int beam_id,int* beam_node_ids, int beam_matice_id, int beam_matpar_id, int beam_numpar_id, ElementProperties* beam_properties);
    43                 Beam(int beam_id,Hook* beam_hnodes, Hook* beam_hmatice, Hook* beam_hmatpar, Hook* beam_hnumpar, ElementProperties* beam_properties);
     39                Beam(int beam_id,int* beam_node_ids, int beam_matice_id, int beam_matpar_id, int beam_numpar_id);
     40                Beam(int beam_id,Hook* beam_hnodes, Hook* beam_hmatice, Hook* beam_hmatpar, Hook* beam_hnumpar, Inputs* beam_inputs);
    4441                Beam(int i, IoModel* iomodel);
    4542                ~Beam();
    … …  
    6158                /*}}}*/
    6259                /*numerics: {{{1*/
    63                 void  CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    64                 void  CreatePVector(Vec pg, void* inputs, int analysis_type,int sub_analysis_type);
    65                 void  UpdateFromInputs(void* inputs);
     60                void  CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type);
     61                void  CreatePVector(Vec pg,  int analysis_type,int sub_analysis_type);
    6662                void  UpdateInputs(double* solution, int analysis_type, int sub_analysis_type);
    6763                void  GetDofList(int* doflist,int* pnumberofdofs);
    6864                void  GetDofList1(int* doflist);
    69                 void  CreateKMatrixDiagnosticHutter(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
     65                void  CreateKMatrixDiagnosticHutter(Mat Kgg,int analysis_type,int sub_analysis_type);
    7066                void  GetParameterValue(double* pp, double* plist, double* gauss_l1l2l3);
    71                 void  CreatePVectorDiagnosticHutter(Vec pg,void* inputs,int analysis_type,int sub_analysis_type);
     67                void  CreatePVectorDiagnosticHutter(Vec pg,int analysis_type,int sub_analysis_type);
    7268                void* GetMatPar();
    7369
    74                 void  ComputeBasalStress(Vec sigma_b,void* inputs,int analysis_type,int sub_analysis_type);
    75                 void  ComputePressure(Vec p_g,void* inputs,int analysis_type,int sub_analysis_type);
    76                 void  ComputeStrainRate(Vec eps,void* inputs,int analysis_type,int sub_analysis_type);
     70                void  ComputeBasalStress(Vec sigma_b,int analysis_type,int sub_analysis_type);
     71                void  ComputePressure(Vec p_g,int analysis_type,int sub_analysis_type);
     72                void  ComputeStrainRate(Vec eps,int analysis_type,int sub_analysis_type);
    7773                void  GetNodes(void** vpnodes);
    7874                /*}}}*/
    7975                /*not implemented: {{{1*/
    80                 int   GetShelf();
    81                 int   GetOnBed();
     76                bool  GetShelf();
     77                bool  GetOnBed();
    8278                void  GetBedList(double*);
    8379                void  GetThicknessList(double* thickness_list);

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

    r3599 r3612  
    1717               
    1818                virtual        ~Element(){};
     19                virtual int    Enum()=0;
    1920                virtual void   Configure(void* loads,void* nodes,void* materials,void* parameters)=0;
    20                 virtual void   CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type)=0;
    21                 virtual void   CreatePVector(Vec pg, void* inputs, int analysis_type,int sub_analysis_type)=0;
    22                 virtual int    Enum()=0;
    23                 virtual void   UpdateFromInputs(void* inputs)=0;
     21               
     22                virtual void   CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type)=0;
     23                virtual void   CreatePVector(Vec pg,  int analysis_type,int sub_analysis_type)=0;
    2424                virtual void   UpdateInputs(double* solution, int analysis_type, int sub_analysis_type)=0;
    2525                virtual void   GetNodes(void** nodes)=0;
    2626                virtual void*  GetMatPar()=0;
    27                 virtual int    GetShelf()=0;
    28                 virtual int    GetOnBed()=0;
     27                virtual bool   GetShelf()=0;
     28                virtual bool   GetOnBed()=0;
    2929                virtual void   GetThicknessList(double* thickness_list)=0;
    3030                virtual void   GetBedList(double* bed_list)=0;
    31                 virtual void   Du(Vec du_g,void* inputs,int analysis_type,int sub_analysis_type)=0;
    32                 virtual void   Gradj(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type,char* control_type)=0;
    33                 virtual void   GradjDrag(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type)=0;
    34                 virtual void   GradjB(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type)=0;
    35                 virtual double Misfit(void* inputs,int analysis_type,int sub_analysis_type)=0;
    36                 virtual double CostFunction(void* inputs,int analysis_type,int sub_analysis_type)=0;
    37                 virtual double SurfaceArea(void* inputs,int analysis_type,int sub_analysis_type)=0;
    38                 virtual void   ComputeBasalStress(Vec sigma_b,void* inputs,int analysis_type,int sub_analysis_type)=0;
    39                 virtual void   ComputePressure(Vec p_g,       void* inputs,int analysis_type,int sub_analysis_type)=0;
    40                 virtual void   ComputeStrainRate(Vec eps,     void* inputs,int analysis_type,int sub_analysis_type)=0;
     31                virtual void   Du(Vec du_g,int analysis_type,int sub_analysis_type)=0;
     32                virtual void   Gradj(Vec grad_g,int analysis_type,int sub_analysis_type,char* control_type)=0;
     33                virtual void   GradjDrag(Vec grad_g,int analysis_type,int sub_analysis_type)=0;
     34                virtual void   GradjB(Vec grad_g,int analysis_type,int sub_analysis_type)=0;
     35                virtual double Misfit(int analysis_type,int sub_analysis_type)=0;
     36                virtual double CostFunction(int analysis_type,int sub_analysis_type)=0;
     37                virtual double SurfaceArea(int analysis_type,int sub_analysis_type)=0;
     38                virtual void   ComputeBasalStress(Vec sigma_b,int analysis_type,int sub_analysis_type)=0;
     39                virtual void   ComputePressure(Vec p_g,       int analysis_type,int sub_analysis_type)=0;
     40                virtual void   ComputeStrainRate(Vec eps,     int analysis_type,int sub_analysis_type)=0;
    4141                virtual double MassFlux(double* segment,double* ug)=0;
    4242

  • issm/trunk/src/c/objects/Friction.cpp

    r2907 r3612  
    4545        printf("\n");
    4646
    47         printf("velocities: ");
    48         for(i=0;i<3;i++)printf("%g ",friction->velocities[i*2+0]);
     47        printf("vx: ");
     48        for(i=0;i<3;i++)printf("%g ",friction->vx[i]);
    4949        printf("\n            ");
    50         for(i=0;i<3;i++)printf("%g ",friction->velocities[i*2+1]);
     50        printf("vy: ");
     51        for(i=0;i<3;i++)printf("%g ",friction->vy[i]);
    5152        printf("\n");
     53        if(friction->vz){
     54                printf("vz: ");
     55                for(i=0;i<3;i++)printf("%g ",friction->vz[i]);
     56                printf("\n");
     57        }
     58
    5259}
    5360/*}}}*/
    … …  
    6673        friction->bed=NULL;
    6774        friction->thickness=NULL;
    68         friction->velocities=NULL;
     75        friction->vx=NULL;
     76        friction->vy=NULL;
     77        friction->vz=NULL;
    6978        friction->p=UNDEF;
    7079        friction->q=UNDEF;
    … …  
    127136                //We need the velocity magnitude to evaluate the basal stress:
    128137                if(strcmp(friction->element_type,"2d")){
    129                         velocity_x[i]=*(friction->velocities+2*i+0); //velocities of size numgridsx2
    130                         velocity_y[i]=*(friction->velocities+2*i+1);
     138                        velocity_x[i]=friction->vx[i];//velocities of size numgridsx2
     139                        velocity_y[i]=friction->vy[i];
    131140                        velocity_mag[i]=sqrt(pow(velocity_x[i],2)+pow(velocity_y[i],2));
    132141                }
    133142                else{
    134                         velocity_x[i]=*(friction->velocities+3*i+0); //velocities of size numgridsx3
    135                         velocity_y[i]=*(friction->velocities+3*i+1);
    136                         velocity_z[i]=*(friction->velocities+3*i+2);
     143                        velocity_x[i]=friction->vx[i]; //velocities of size numgridsx3
     144                        velocity_y[i]=friction->vy[i];
     145                        velocity_z[i]=friction->vz[i];
    137146                        velocity_mag[i]=sqrt(pow(velocity_x[i],2)+pow(velocity_y[i],2)+pow(velocity_z[i],2));
    138147                }
    … …  
    177186
    178187                //We need the velocity magnitude to evaluate the basal stress:
    179                 velocity_x[i]=*(friction->velocities+2*i+0); //velocities of size numgridsx2
    180                 velocity_y[i]=*(friction->velocities+2*i+1);
     188                velocity_x[i]=friction->vx[i]; //velocities of size numgridsx2
     189                velocity_y[i]=friction->vy[i];
    181190                velocity_mag[i]=sqrt(pow(velocity_x[i],2)+pow(velocity_y[i],2));
    182191       

  • issm/trunk/src/c/objects/Friction.h

    r2359 r3612  
    1616        double* bed;
    1717        double* thickness;
    18         double* velocities;
     18        double* vx;
     19        double* vy;
     20        double* vz;
    1921        double  p;
    2022        double  q;

  • issm/trunk/src/c/objects/Hook.cpp

    r3570 r3612  
    1616#include "./Hook.h"
    1717#include "../EnumDefinitions/EnumDefinitions.h"
    18 #include "./ParameterInputs.h"
    1918#include "../shared/shared.h"
    2019#include "../include/typedefs.h"

  • issm/trunk/src/c/objects/Input.h

    r3463 r3612  
    1 /*!\file Input.h
    2  * \brief: header file for input object
    3  */
     1/*!\file: Input.h
     2 * \brief abstract class for Input object
     3 */ 
    44
    5 #ifndef _INPUT_H_
    6 #define _INPUT_H_
    75
    8 #include "./Param.h" //for enums STRING,INTEGER,DOUBLE AND DOUBLEVEC
    9 #include "./Node.h"
     6#ifndef _EINPUT_H_
     7#define _EINPUT_H_
    108
    11 #define INPUTSTRING 200 //max string length
     9#include "./Object.h"
    1210
    13 class Input: public Object {
     11class Input: public Object{
    1412
    15         private:
    16                 char name[INPUTSTRING]; //"analysis_type","velocity", etc ...
    17                 int  type; //STRING, INTEGER, DOUBLE OR DOUBLEVEC
    18 
    19                 /*placeholders: */
    20                 char string[INPUTSTRING];
    21                 int  integer;
    22                 double scalar;
    23                 double* vector;
     13        public:
    2414               
    25                 int  ndof; //number of dofs for input array
    26                 int  numberofnodes; //size of input array for each dof
    27 
    28         public:
    29 
    30                 /*constructors and destructors: */
    31                 Input();
    32                 ~Input();
    33                 Input(char* name,char* string);
    34                 Input(char* name,int integer);
    35                 Input(char* name,double scalar);
    36                 Input(char* name,double* vector,int ndof,int numberofnodes);
    37                 Input(char* name,Vec petscvector,int ndof, int numberofnodes);
    38 
    39                
    40                 Object* copy();
    41 
    42                 /*fill virtual routines: */
    43                 void  Echo();
    44                 void  DeepEcho(void);
    45                 void  Marshall(char** pmarshalled_dataset);
    46                 int   MarshallSize();
    47                 char* GetName();
    48                 void  Demarshall(char** pmarshalled_dataset);
    49                 int   Enum();
    50                 int   GetId();
    51                 int   MyRank();
    52 
    53                 /*Input specific routines: */
    54                 int   IsSameName(char* input_name);
    55                 int   IsPresent(char* input_name);
    56                 void  Recover(double* values, int ndof, int* dofs,int numnodes,void** nodes);
    57                 void  Recover(int* pinteger);
    58                 void  Recover(char** pstring);
    59                 void  Recover(double* pdouble);
    60 
    61                 /*get into workspace: */
    62                 Vec   Get(int* dofs, int numdofs);
     15                virtual        ~Input(){};
     16                virtual int    Enum()=0; //object Enum
     17                virtual int    EnumType()=0; //type of input (vx,vy?)
    6318
    6419};
    65 
    66 #endif  /* _INPUT_H_ */
     20#endif

  • issm/trunk/src/c/objects/Matice.cpp

    r3595 r3612  
    5353        /*Compute B on the element if provided*/
    5454        if (num_vertices==3 || num_vertices==6){
    55                 if (iomodel->B){
     55                if (iomodel->rheology_B){
    5656                        B_avg=0;
    5757                        for(j=0;j<num_vertices;j++){
    58                                 B_avg+=*(iomodel->B+((int)*(iomodel->elements+num_vertices*i+j)-1));
     58                                B_avg+=*(iomodel->rheology_B+((int)*(iomodel->elements+num_vertices*i+j)-1));
    5959                        }
    6060                        B_avg=B_avg/num_vertices;
    … …  
    6262        }
    6363        else if (num_vertices==1 || num_vertices==2){
    64                 if (iomodel->B){
    65                         B_avg=*(iomodel->B+i);
     64                if (iomodel->rheology_B){
     65                        B_avg=*(iomodel->rheology_B+i);
    6666                }
    6767        }
    6868        else ISSMERROR("num_vertices = %i not supported yet",num_vertices);
    6969       
    70         if (iomodel->B) matice_B=B_avg;
     70        if (iomodel->rheology_B) matice_B=B_avg;
    7171        else            matice_B=UNDEF;
    72         if (iomodel->n){
     72        if (iomodel->rheology_n){
    7373                if (num_vertices==3 || num_vertices==6){
    74                         matice_n=(double)*(iomodel->n+i);
     74                        matice_n=(double)*(iomodel->rheology_n+i);
    7575                }
    7676                else if (num_vertices==1 || num_vertices==2){
    7777                        /*n is on the elements for now, so just take the first one*/
    78                         matice_n=(double)*(iomodel->n);
     78                        matice_n=(double)*(iomodel->rheology_n);
    7979                }
    8080                else ISSMERROR("num_vertices = %i not supported yet",num_vertices);

  • issm/trunk/src/c/objects/Node.cpp

    r3570 r3612  
    1717#include <string.h>
    1818#include "../EnumDefinitions/EnumDefinitions.h"
    19 #include "./ParameterInputs.h"
    2019#include "../shared/shared.h"
    2120#include "../include/typedefs.h"

  • issm/trunk/src/c/objects/Numpar.cpp

    r3607 r3612  
    3939        cm_maxdmp_value=UNDEF;
    4040        cm_maxdmp_slope=UNDEF;
     41        dt=UNDEF;
    4142
    4243        return;
    … …  
    8283        memcpy(marshalled_dataset,&cm_maxdmp_value,sizeof(cm_maxdmp_value));marshalled_dataset+=sizeof(cm_maxdmp_value);
    8384        memcpy(marshalled_dataset,&cm_maxdmp_slope,sizeof(cm_maxdmp_slope));marshalled_dataset+=sizeof(cm_maxdmp_slope);
     85        memcpy(marshalled_dataset,&dt,sizeof(dt));marshalled_dataset+=sizeof(dt);
    8486       
    8587        *pmarshalled_dataset=marshalled_dataset;
    … …  
    101103                +sizeof(cm_maxdmp_value)
    102104                +sizeof(cm_maxdmp_slope)
     105                +sizeof(dt)
    103106                +sizeof(int); //sizeof(int) for enum type
    104107}
    … …  
    128131        memcpy(&cm_maxdmp_value,marshalled_dataset,sizeof(cm_maxdmp_value));marshalled_dataset+=sizeof(cm_maxdmp_value);
    129132        memcpy(&cm_maxdmp_slope,marshalled_dataset,sizeof(cm_maxdmp_slope));marshalled_dataset+=sizeof(cm_maxdmp_slope);
     133        memcpy(&dt,marshalled_dataset,sizeof(dt));marshalled_dataset+=sizeof(dt);
    130134
    131135        /*return: */
    … …  
    158162        parameters->FindParam(&cm_maxdmp_value,"cm_maxdmp_value");
    159163        parameters->FindParam(&cm_maxdmp_slope,"cm_maxdmp_slope");
     164        parameters->FindParam(&dt,"dt");
    160165
    161166        return;
    … …  
    185190        printf("   cm_maxdmp_value: %g\n",cm_maxdmp_value);
    186191        printf("   cm_maxdmp_slope: %g\n",cm_maxdmp_slope);
     192        printf("   dt: %g\n",dt);
    187193}
    188194/*}}}*/
    … …  
    203209        printf("   cm_maxdmp_value: %g\n",cm_maxdmp_value);
    204210        printf("   cm_maxdmp_slope: %g\n",cm_maxdmp_slope);
     211        printf("   dt: %g\n",dt);
    205212}
    206213/*}}}*/
    … …  
    246253        inputs->Recover("cm_maxdmp_value",&cm_maxdmp_value);
    247254        inputs->Recover("cm_maxdmp_slope",&cm_maxdmp_slope);
    248 
    249 }
    250 /*}}}*/
     255        inputs->Recover("dt",&dt);
     256
     257}
     258/*}}}*/

  • issm/trunk/src/c/objects/Numpar.h

    r3463 r3612  
    2727                double cm_maxdmp_value;
    2828                double cm_maxdmp_slope;
     29                double dt;
    2930
    3031                Numpar();

  • issm/trunk/src/c/objects/OptArgs.h

    r1881 r3612  
    2424#else
    2525
    26 #include "./ParameterInputs.h"
    27 class ParameterInputs;
    28 
    2926struct OptArgs{
    3027        Model* model;
    3128        double* param_g;
    3229        double* grad_g;
    33         ParameterInputs* inputs;
    3430        int n;
    3531};

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

    r3599 r3612  
    266266}
    267267/*}}}*/
    268 /*FUNCTION UpdateFromInputs {{{1*/
    269 void  Penta::UpdateFromInputs(void* vinputs){
    270 
    271         int     dofs[1]={0};
    272         double  temperature_list[6];
    273         double  temperature_average;
    274         double  B_list[6];
    275         double  B_average;
    276 
    277         /*dynamic objects pointed to by hooks: */
    278         Node**  nodes=NULL;
    279         Matice* matice=NULL;
    280 
    281         ParameterInputs* inputs=NULL;
    282 
    283         /*If on water, skip: */
    284         if(this->properties.onwater)return;
    285 
    286         /*recover objects from hooks: */
    287         nodes=(Node**)hnodes.deliverp();
    288         matice=(Matice*)hmatice.delivers();
    289 
    290         /*recover pointers: */
    291         inputs=(ParameterInputs*)vinputs;
    292 
    293         /*Update internal data if inputs holds new values: */
    294         if (this->properties.h)              inputs->Recover("thickness",&this->properties.h[0],1,dofs,6,(void**)nodes);
    295         if (this->properties.s)              inputs->Recover("surface",&this->properties.s[0],1,dofs,6,(void**)nodes);
    296         if (this->properties.b)              inputs->Recover("bed",&this->properties.b[0],1,dofs,6,(void**)nodes);
    297         if (this->properties.k)              inputs->Recover("drag",&this->properties.k[0],1,dofs,6,(void**)nodes);
    298         if (this->properties.melting)        inputs->Recover("melting",&this->properties.melting[0],1,dofs,6,(void**)nodes);
    299         if (this->properties.accumulation)   inputs->Recover("accumulation",&this->properties.accumulation[0],1,dofs,6,(void**)nodes);
    300         if (this->properties.geothermalflux) inputs->Recover("geothermalflux",&this->properties.geothermalflux[0],1,dofs,6,(void**)nodes);
    301 
    302         //Update material if necessary
    303         if(inputs->Recover("temperature",&temperature_list[0],1,dofs,6,(void**)nodes)){
    304                 if(matice && !this->properties.collapse){
    305                         //B_average=(Paterson(temperature_list[0])+Paterson(temperature_list[1])+Paterson(temperature_list[2])
    306                         //                      +Paterson(temperature_list[3])+Paterson(temperature_list[4])+Paterson(temperature_list[5]))/6.0;
    307                         temperature_average=(temperature_list[0]+temperature_list[1]+temperature_list[2]+temperature_list[3]+temperature_list[4]+temperature_list[5])/6.0;
    308                         B_average=Paterson(temperature_average);
    309                         matice->SetB(B_average);
    310                 }
    311         }
    312 
    313         if(inputs->Recover("temperature_average",&temperature_list[0],1,dofs,6,(void**)nodes)){
    314                 if(matice && this->properties.collapse){
    315                         temperature_average=(temperature_list[0]+temperature_list[1]+temperature_list[2]+temperature_list[3]+temperature_list[4]+temperature_list[5])/6.0;
    316                         B_average=Paterson(temperature_average);
    317                         //B_average=(Paterson(temperature_list[0])+Paterson(temperature_list[1])+Paterson(temperature_list[2])
    318                         //                      +Paterson(temperature_list[3])+Paterson(temperature_list[4])+Paterson(temperature_list[5]))/6.0;
    319                         matice->SetB(B_average);
    320                 }
    321         }
    322 
    323         if(inputs->Recover("B",&B_list[0],1,dofs,6,(void**)nodes)){
    324                 if(matice){
    325                         B_average=(B_list[0]+B_list[1]+B_list[2]+B_list[3]+B_list[4]+B_list[5])/6.0;
    326                         matice->SetB(B_average);
    327                 }
    328         }
    329 
    330 }
    331 /*}}}*/
    332268/*FUNCTION UpdateFromDakota {{{1*/
    333269void  Penta::UpdateFromDakota(void* vinputs){
    … …  
    345281/*Object functions*/
    346282/*FUNCTION ComputeBasalStress {{{1*/
    347 void  Penta::ComputeBasalStress(Vec sigma_b,void* vinputs,int analysis_type,int sub_analysis_type){
     283void  Penta::ComputeBasalStress(Vec sigma_b,int analysis_type,int sub_analysis_type){
    348284
    349285        int i,j;
    … …  
    477413/*}}}*/
    478414/*FUNCTION ComputePressure {{{1*/
    479 void  Penta::ComputePressure(Vec pg,void* vinputs,int analysis_type,int sub_analysis_type){
     415void  Penta::ComputePressure(Vec pg,int analysis_type,int sub_analysis_type){
    480416
    481417        int i;
    … …  
    519455/*}}}*/
    520456/*FUNCTION ComputeStrainRate {{{1*/
    521 void  Penta::ComputeStrainRate(Vec eps,void* vinputs,int analysis_type,int sub_analysis_type){
     457void  Penta::ComputeStrainRate(Vec eps,int analysis_type,int sub_analysis_type){
    522458
    523459        ISSMERROR("Not implemented yet");
    … …  
    526462/*}}}*/
    527463/*FUNCTION CostFunction {{{1*/
    528 double Penta::CostFunction(void* inputs,int analysis_type,int sub_analysis_type){
     464double Penta::CostFunction(int analysis_type,int sub_analysis_type){
    529465
    530466        double J;
    … …  
    560496/*FUNCTION CreateKMatrix {{{1*/
    561497
    562 void  Penta::CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     498void  Penta::CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type){
    563499
    564500        /*Just branch to the correct element stiffness matrix generator, according to the type of analysis we are carrying out: */
    … …  
    608544/*}}}*/
    609545/*FUNCTION CreateKMatrixDiagnosticHoriz {{{1*/
    610 void Penta::CreateKMatrixDiagnosticHoriz( Mat Kgg, void* vinputs, int analysis_type,int sub_analysis_type){
     546void Penta::CreateKMatrixDiagnosticHoriz( Mat Kgg, int analysis_type,int sub_analysis_type){
    611547
    612548
    … …  
    844780/*}}}*/
    845781/*FUNCTION CreateKMatrixDiagnosticStokes {{{1*/
    846 void Penta::CreateKMatrixDiagnosticStokes( Mat Kgg, void* vinputs, int analysis_type,int sub_analysis_type){
     782void Penta::CreateKMatrixDiagnosticStokes( Mat Kgg, int analysis_type,int sub_analysis_type){
    847783
    848784        int i,j;
    … …  
    11391075/*}}}*/
    11401076/*FUNCTION CreateKMatrixDiagnosticVert {{{1*/
    1141 void Penta::CreateKMatrixDiagnosticVert( Mat Kgg, void* vinputs, int analysis_type,int sub_analysis_type){
     1077void Penta::CreateKMatrixDiagnosticVert( Mat Kgg, int analysis_type,int sub_analysis_type){
    11421078
    11431079        /* local declarations */
    … …  
    12751211/*}}}*/
    12761212/*FUNCTION CreateKMatrixMelting {{{1*/
    1277 void  Penta::CreateKMatrixMelting(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     1213void  Penta::CreateKMatrixMelting(Mat Kgg,int analysis_type,int sub_analysis_type){
    12781214
    12791215        Tria* tria=NULL;
    … …  
    12961232/*FUNCTION CreateKMatrixPrognostic {{{1*/
    12971233
    1298 void  Penta::CreateKMatrixPrognostic(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     1234void  Penta::CreateKMatrixPrognostic(Mat Kgg,int analysis_type,int sub_analysis_type){
    12991235
    13001236        /*Collapsed formulation: */
    … …  
    13171253/*FUNCTION CreateKMatrixSlopeCompute {{{1*/
    13181254
    1319 void  Penta::CreateKMatrixSlopeCompute(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     1255void  Penta::CreateKMatrixSlopeCompute(Mat Kgg,int analysis_type,int sub_analysis_type){
    13201256
    13211257        /*Collapsed formulation: */
    … …  
    13371273/*}}}*/
    13381274/*FUNCTION CreateKMatrixThermal {{{1*/
    1339 void  Penta::CreateKMatrixThermal(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     1275void  Penta::CreateKMatrixThermal(Mat Kgg,int analysis_type,int sub_analysis_type){
    13401276
    13411277        /* local declarations */
    … …  
    16001536/*}}}*/
    16011537/*FUNCTION CreatePVector {{{1*/
    1602 void  Penta::CreatePVector(Vec pg, void* inputs, int analysis_type,int sub_analysis_type){
     1538void  Penta::CreatePVector(Vec pg, int analysis_type,int sub_analysis_type){
    16031539
    16041540        /*Just branch to the correct element stiffness matrix generator, according to the type of analysis we are carrying out: */
    … …  
    16461582/*}}}*/
    16471583/*FUNCTION CreatePVectorDiagnosticHoriz {{{1*/
    1648 void Penta::CreatePVectorDiagnosticHoriz( Vec pg, void* vinputs,int analysis_type,int sub_analysis_type){
     1584void Penta::CreatePVectorDiagnosticHoriz( Vec pg, int analysis_type,int sub_analysis_type){
    16491585
    16501586        int i,j;
    … …  
    18011737/*}}}*/
    18021738/*FUNCTION CreatePVectorDiagnosticStokes {{{1*/
    1803 void Penta::CreatePVectorDiagnosticStokes( Vec pg, void* vinputs,int analysis_type,int sub_analysis_type){
     1739void Penta::CreatePVectorDiagnosticStokes( Vec pg, int analysis_type,int sub_analysis_type){
    18041740
    18051741        /*indexing: */
    … …  
    20511987/*}}}*/
    20521988/*FUNCTION CreatePVectorDiagnosticVert {{{1*/
    2053 void  Penta::CreatePVectorDiagnosticVert( Vec pg, void* vinputs,int analysis_type,int sub_analysis_type){
     1989void  Penta::CreatePVectorDiagnosticVert( Vec pg, int analysis_type,int sub_analysis_type){
    20541990
    20551991        int i;
    … …  
    21902126/*}}}*/
    21912127/*FUNCTION CreatePVectorMelting {{{1*/
    2192 void Penta::CreatePVectorMelting( Vec pg, void* vinputs,int analysis_type,int sub_analysis_type){
     2128void Penta::CreatePVectorMelting( Vec pg, int analysis_type,int sub_analysis_type){
    21932129        return;
    21942130}
    … …  
    21962132/*FUNCTION CreatePVectorPrognostic {{{1*/
    21972133
    2198 void Penta::CreatePVectorPrognostic( Vec pg, void* inputs, int analysis_type,int sub_analysis_type){
     2134void Penta::CreatePVectorPrognostic( Vec pg, int analysis_type,int sub_analysis_type){
    21992135
    22002136        /*Collapsed formulation: */
    … …  
    22162152/*FUNCTION CreatePVectorSlopeCompute {{{1*/
    22172153
    2218 void Penta::CreatePVectorSlopeCompute( Vec pg, void* inputs, int analysis_type,int sub_analysis_type){
     2154void Penta::CreatePVectorSlopeCompute( Vec pg, int analysis_type,int sub_analysis_type){
    22192155
    22202156        /*Collapsed formulation: */
    … …  
    22352171/*}}}*/
    22362172/*FUNCTION CreatePVectorThermal {{{1*/
    2237 void Penta::CreatePVectorThermal( Vec pg, void* vinputs,int analysis_type,int sub_analysis_type){
     2173void Penta::CreatePVectorThermal( Vec pg, int analysis_type,int sub_analysis_type){
    22382174
    22392175
    … …  
    24382374/*}}}*/
    24392375/*FUNCTION Du {{{1*/
    2440 void  Penta::Du(Vec du_g,void* inputs,int analysis_type,int sub_analysis_type){
     2376void  Penta::Du(Vec du_g,int analysis_type,int sub_analysis_type){
    24412377
    24422378        int i;
    … …  
    36473583/*}}}*/
    36483584/*FUNCTION GetOnBed {{{1*/
    3649 int Penta::GetOnBed(){
     3585bool Penta::GetOnBed(){
    36503586        return this->properties.onbed;
    36513587}
    … …  
    37923728/*}}}*/
    37933729/*FUNCTION Gradj {{{1*/
    3794 void  Penta::Gradj(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type,char* control_type){
     3730void  Penta::Gradj(Vec grad_g,int analysis_type,int sub_analysis_type,char* control_type){
    37953731
    37963732        /*If on water, skip grad (=0): */
    … …  
    38073743/*}}}*/
    38083744/*FUNCTION GradjDrag {{{1*/
    3809 void  Penta::GradjDrag(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type){
     3745void  Penta::GradjDrag(Vec grad_g,int analysis_type,int sub_analysis_type){
    38103746
    38113747        Tria* tria=NULL;
    … …  
    38403776/*}}}*/
    38413777/*FUNCTION GradjB {{{1*/
    3842 void  Penta::GradjB(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type){
     3778void  Penta::GradjB(Vec grad_g,int analysis_type,int sub_analysis_type){
    38433779
    38443780        Tria* tria=NULL;
    … …  
    38733809/*}}}*/
    38743810/*FUNCTION Misfit {{{1*/
    3875 double Penta::Misfit(void* inputs,int analysis_type,int sub_analysis_type){
     3811double Penta::Misfit(int analysis_type,int sub_analysis_type){
    38763812
    38773813        double J;
    … …  
    40063942/*}}}1*/
    40073943/*FUNCTION SurfaceArea {{{1*/
    4008 double Penta::SurfaceArea(void* inputs,int analysis_type,int sub_analysis_type){
     3944double Penta::SurfaceArea(int analysis_type,int sub_analysis_type){
    40093945
    40103946        double S;

  • issm/trunk/src/c/objects/Penta.h

    r3599 r3612  
    99class Node;
    1010class Hook;
    11 class ElementProperties;
    1211class DataSet;
    1312struct IoModel;
    … …  
    2120#include "./Tria.h"
    2221#include "./Hook.h"
    23 #include "./ElementProperties.h"
    2422#include "../ModelProcessorx/IoModel.h"
    25 #include "./ParameterInputs.h"
    2623#include "./Node.h"
    2724
    … …  
    3633                Hook hnumpar; //hook to 1 numpar
    3734
    38                 ElementProperties properties;
    39                 DataSet* inputs;
     35                Inputs* inputs;
    4036
    4137        public:
    … …  
    4339                /*FUNCTION constructors, destructors {{{1*/
    4440                Penta();
    45                 Penta(int penta_id,int* penta_node_ids, int penta_matice_id, int penta_matpar_id, int penta_numpar_id, ElementProperties* properties);
    46                 Penta(int penta_id,Hook* penta_hnodes, Hook* penta_hmatice, Hook* penta_hmatpar, Hook* penta_hnumpar, ElementProperties* penta_properties);
     41                Penta(int penta_id,int* penta_node_ids, int penta_matice_id, int penta_matpar_id, int penta_numpar_id);
     42                Penta(int penta_id,Hook* penta_hnodes, Hook* penta_hmatice, Hook* penta_hmatpar, Hook* penta_hnumpar, Inputs* inputs);
    4743                Penta(int i, IoModel* iomodel);
    4844                ~Penta();
    … …  
    6258                void*  SpawnTria(int g0, int g1, int g2);
    6359                void  UpdateFromDakota(void* inputs);
    64                 void  UpdateFromInputs(void* inputs);
    6560                void  UpdateInputs(double* solution, int analysis_type, int sub_analysis_type);
    6661                void  SetClone(int* minranks);
    … …  
    6863                /*}}}*/
    6964                /*FUNCTION element numerical routines {{{1*/
    70                 void  CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    71                 void  CreateKMatrixDiagnosticHoriz( Mat Kgg, void* inputs, int analysis_type,int sub_analysis_type);
    72                 void  CreateKMatrixDiagnosticVert( Mat Kgg, void* inputs, int analysis_type,int sub_analysis_type);
    73                 void  CreatePVector(Vec pg, void* inputs, int analysis_type,int sub_analysis_type);
     65                void  CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type);
     66                void  CreateKMatrixDiagnosticHoriz( Mat Kgg, int analysis_type,int sub_analysis_type);
     67                void  CreateKMatrixDiagnosticVert( Mat Kgg, int analysis_type,int sub_analysis_type);
     68                void  CreatePVector(Vec pg, int analysis_type,int sub_analysis_type);
    7469                void  GetDofList(int* doflist,int* pnumberofdofs);
    7570                void  GetDofList1(int* doflist);
    7671                void* GetMatPar();
    77                 int   GetShelf();
     72                bool   GetShelf();
    7873                void  GetNodes(void** nodes);
    79                 int   GetOnBed();
    80                 void  Du(Vec du_g,void* inputs,int analysis_type,int sub_analysis_type);
    81                 void  Gradj(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type,char* control_type);
    82                 void  GradjDrag(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type);
    83                 void  GradjB(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type);
    84                 double Misfit(void* inputs,int analysis_type,int sub_analysis_type);
    85                 double SurfaceArea(void* inputs,int analysis_type,int sub_analysis_type);
    86                 double CostFunction(void* inputs,int analysis_type,int sub_analysis_type);
     74                bool   GetOnBed();
     75                void  Du(Vec du_g,int analysis_type,int sub_analysis_type);
     76                void  Gradj(Vec grad_g,int analysis_type,int sub_analysis_type,char* control_type);
     77                void  GradjDrag(Vec grad_g,int analysis_type,int sub_analysis_type);
     78                void  GradjB(Vec grad_g,int analysis_type,int sub_analysis_type);
     79                double Misfit(int analysis_type,int sub_analysis_type);
     80                double SurfaceArea(int analysis_type,int sub_analysis_type);
     81                double CostFunction(int analysis_type,int sub_analysis_type);
    8782               
    8883                void          GetThicknessList(double* thickness_list);
    … …  
    9994                void  GetNodalFunctionsDerivativesReference(double* dl1dl6,double* gauss_coord);
    10095                void  GetJacobianInvert(double*  Jinv, double* xyz_list,double* gauss_coord);
    101                 void  CreatePVectorDiagnosticHoriz( Vec pg, void* inputs,int analysis_type,int sub_analysis_type);
    102                 void  CreatePVectorDiagnosticVert( Vec pg, void* inputs,int analysis_type,int sub_analysis_type);
     96                void  CreatePVectorDiagnosticHoriz( Vec pg, int analysis_type,int sub_analysis_type);
     97                void  CreatePVectorDiagnosticVert( Vec pg, int analysis_type,int sub_analysis_type);
    10398                void  GetParameterValue(double* pvalue, double* v_list,double* gauss_coord);
    10499                void  GetParameterDerivativeValue(double* p, double* p_list,double* xyz_list, double* gauss_coord);
    105100                void  GetNodalFunctions(double* l1l6, double* gauss_coord);
    106101                void  FieldExtrude(Vec field,double* field_serial,char* field_name, int iscollapsed);
    107                 void  ComputeBasalStress(Vec sigma_b,void* inputs,int analysis_type,int sub_analysis_type);
    108                 void  ComputePressure(Vec p_g,void* inputs,int analysis_type,int sub_analysis_type);
    109                 void  ComputeStrainRate(Vec eps,void* inputs,int analysis_type,int sub_analysis_type);
    110                 void  CreateKMatrixSlopeCompute(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type);
    111                 void  CreatePVectorSlopeCompute( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type);
    112                 void  CreateKMatrixPrognostic(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type);
    113                 void  CreatePVectorPrognostic( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type);
     102                void  ComputeBasalStress(Vec sigma_b,int analysis_type,int sub_analysis_type);
     103                void  ComputePressure(Vec p_g,int analysis_type,int sub_analysis_type);
     104                void  ComputeStrainRate(Vec eps,int analysis_type,int sub_analysis_type);
     105                void  CreateKMatrixSlopeCompute(Mat Kgg,int analysis_type,int sub_analysis_type);
     106                void  CreatePVectorSlopeCompute( Vec pg, int analysis_type,int sub_analysis_type);
     107                void  CreateKMatrixPrognostic(Mat Kgg,int analysis_type,int sub_analysis_type);
     108                void  CreatePVectorPrognostic( Vec pg, int analysis_type,int sub_analysis_type);
    114109
    115                 void  CreateKMatrixDiagnosticStokes( Mat Kgg, void* vinputs, int analysis_type,int sub_analysis_type);
    116                 void  CreatePVectorDiagnosticStokes( Vec pg, void* vinputs,int analysis_type,int sub_analysis_type);
     110                void  CreateKMatrixDiagnosticStokes( Mat Kgg, int analysis_type,int sub_analysis_type);
     111                void  CreatePVectorDiagnosticStokes( Vec pg, int analysis_type,int sub_analysis_type);
    117112                void  ReduceMatrixStokes(double* Ke_reduced, double* Ke_temp);
    118113                void  GetMatrixInvert(double*  Ke_invert, double* Ke);
    … …  
    127122                void  ReduceVectorStokes(double* Pe_reduced, double* Ke_temp, double* Pe_temp);
    128123                void  GetNodalFunctionsStokes(double* l1l7, double* gauss_coord);
    129                 void  CreateKMatrixThermal(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
     124                void  CreateKMatrixThermal(Mat Kgg,int analysis_type,int sub_analysis_type);
    130125                void  GetB_conduct(double* B_conduct, double* xyz_list, double* gauss_coord);
    131126                void  GetB_advec(double* B_advec, double* xyz_list, double* gauss_coord);
    132127                void  GetBprime_advec(double* Bprime_advec, double* xyz_list, double* gauss_coord);
    133128                void  GetB_artdiff(double* B_artdiff, double* xyz_list, double* gauss_coord);
    134                 void  CreateKMatrixMelting(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    135                 void  CreatePVectorThermal( Vec pg, void* vinputs,int analysis_type,int sub_analysis_type);
    136                 void  CreatePVectorMelting( Vec pg, void* vinputs,int analysis_type,int sub_analysis_type);
     129                void  CreateKMatrixMelting(Mat Kgg,int analysis_type,int sub_analysis_type);
     130                void  CreatePVectorThermal( Vec pg, int analysis_type,int sub_analysis_type);
     131                void  CreatePVectorMelting( Vec pg, int analysis_type,int sub_analysis_type);
    137132                void  GetPhi(double* phi, double*  epsilon, double viscosity);
    138133                double MassFlux(double* segment,double* ug);

  • issm/trunk/src/c/objects/Result.cpp

    r3567 r3612  
    1515#include "../EnumDefinitions/EnumDefinitions.h"
    1616#include "../include/macros.h"
    17 #include "./ParameterInputs.h"
    1817#include "../shared/shared.h"
    1918#include "../include/typedefs.h"

  • issm/trunk/src/c/objects/Riftfront.h

    r3463 r3612  
    1010#include "./Element.h"
    1111#include "./Node.h"
    12 #include "./ParameterInputs.h"
    1312
    1413#define MAX_RIFTFRONT_GRIDS 2 //max number of grids on a rift flank, only 2 because 2d for now.

  • issm/trunk/src/c/objects/Sing.cpp

    r3599 r3612  
    226226}
    227227/*}}}*/
    228 /*FUNCTION Sing::UpdateFromInputs {{{1*/
    229 void  Sing::UpdateFromInputs(void* vinputs){
    230 
    231         int     dofs[1]={0};
    232         double  temperature;
    233         double  B;
    234         double  B_average;
     228/*FUNCTION Sing::UpdateInputs {{{1*/
     229void  Sing::UpdateInputs(double* solution, int analysis_type, int sub_analysis_type){
     230        ISSMERROR(" not supported yet!");
     231}
     232/*}}}*/
     233               
     234/*Object functions*/
     235/*FUNCTION Sing::ComputeBasalStress {{{1*/
     236void  Sing::ComputeBasalStress(Vec p_g,int analysis_type,int sub_analysis_type){
     237
     238        ISSMERROR("Not implemented yet");
     239
     240}
     241/*}}}*/
     242/*FUNCTION Sing::ComputePressure {{{1*/
     243void  Sing::ComputePressure(Vec p_g,int analysis_type,int sub_analysis_type){
     244
     245        int i;
     246        const int numgrids=1;
     247        int doflist[numgrids];
     248        double pressure[numgrids];
     249        double rho_ice,g;
    235250
    236251        /*dynamic objects pointed to by hooks: */
    … …  
    240255        Numpar* numpar=NULL;
    241256
    242         ParameterInputs* inputs=NULL;
    243 
    244         /*recover pointers: */
    245         inputs=(ParameterInputs*)vinputs;
     257        /*Get dof list on which we will plug the pressure values: */
     258        GetDofList1(&doflist[0]);
    246259
    247260        /*recover objects from hooks: */
    … …  
    251264        numpar=(Numpar*)hnumpar.delivers();
    252265
    253         /*Update internal data if inputs holds new values: */
    254         //if (id==1) printf("WARNING if QMU: no hydrostatic equilibrium is applied here (conflict with prognostic, which does not have matpar)\n");
    255         //For now
    256         if(this->properties.h) inputs->Recover("thickness",&this->properties.h[0],1,dofs,1,(void**)nodes);
    257         //Later
    258         /*
    259                 if(inputs->Recover("thickness",&new_h[0],1,dofs,1,(void**)nodes)){
    260         //density, needed later:
    261         double di=(this->matpar->GetRhoIce()/this->matpar->GetRhoWater());
    262         //Go through grids:
    263         for (i=0;i<1;i++){
    264         if(nodes[i]->IsOnShelf()){
    265         this->b[i]=this->b[i]-di*(new_h[i]-h[i]); //hydrostatic equilibrium;
    266         }
    267         this->s[i]=this->b[i]+new_h[i];
    268         this->h[i]=new_h[i];
    269         }
    270         }
    271         */
    272         if (this->properties.k) inputs->Recover("drag",&this->properties.k[0],1,dofs,1,(void**)nodes);
    273        
    274         //Update material if necessary
    275         if(inputs->Recover("temperature_average",&temperature,1,dofs,1,(void**)nodes)){
    276                 B_average=Paterson(temperature);
    277                 matice->SetB(B_average);
    278         }
    279         else if(inputs->Recover("temperature",&temperature,1,dofs,1,(void**)nodes)){
    280                 B=Paterson(temperature);
    281                 matice->SetB(B);
    282         }
    283        
    284         if(inputs->Recover("B",&B,1,dofs,1,(void**)nodes)){
    285                 matice->SetB(B);
    286         }
    287 
    288 }
    289 /*}}}*/
    290 /*FUNCTION Sing::UpdateInputs {{{1*/
    291 void  Sing::UpdateInputs(double* solution, int analysis_type, int sub_analysis_type){
    292         ISSMERROR(" not supported yet!");
    293 }
    294 /*}}}*/
    295                
    296 /*Object functions*/
    297 /*FUNCTION Sing::ComputeBasalStress {{{1*/
    298 void  Sing::ComputeBasalStress(Vec p_g,void* inputs,int analysis_type,int sub_analysis_type){
    299 
    300         ISSMERROR("Not implemented yet");
    301 
    302 }
    303 /*}}}*/
    304 /*FUNCTION Sing::ComputePressure {{{1*/
    305 void  Sing::ComputePressure(Vec p_g,void* inputs,int analysis_type,int sub_analysis_type){
    306 
    307         int i;
    308         const int numgrids=1;
    309         int doflist[numgrids];
    310         double pressure[numgrids];
    311         double rho_ice,g;
    312 
    313         /*dynamic objects pointed to by hooks: */
    314         Node**  nodes=NULL;
    315         Matpar* matpar=NULL;
    316         Matice* matice=NULL;
    317         Numpar* numpar=NULL;
    318 
    319         /*Get dof list on which we will plug the pressure values: */
    320         GetDofList1(&doflist[0]);
    321 
    322         /*recover objects from hooks: */
    323         nodes=(Node**)hnodes.deliverp();
    324         matpar=(Matpar*)hmatpar.delivers();
    325         matice=(Matice*)hmatice.delivers();
    326         numpar=(Numpar*)hnumpar.delivers();
    327 
    328266        /*pressure is lithostatic: */
    329267        rho_ice=matpar->GetRhoIce();
    … …  
    337275/*}}}*/
    338276/*FUNCTION Sing::ComputeStrainRate {{{1*/
    339 void  Sing::ComputeStrainRate(Vec p_g,void* inputs,int analysis_type,int sub_analysis_type){
     277void  Sing::ComputeStrainRate(Vec p_g,int analysis_type,int sub_analysis_type){
    340278
    341279        ISSMERROR("Not implemented yet");
    … …  
    350288/*FUNCTION Sing::CreateKMatrix {{{1*/
    351289
    352 void  Sing::CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     290void  Sing::CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type){
    353291
    354292        /*Just branch to the correct element stiffness matrix generator, according to the type of analysis we are carrying out: */
    … …  
    366304/*FUNCTION Sing::CreateKMatrixDiagnosticHutter {{{1*/
    367305
    368 void  Sing::CreateKMatrixDiagnosticHutter(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     306void  Sing::CreateKMatrixDiagnosticHutter(Mat Kgg,int analysis_type,int sub_analysis_type){
    369307       
    370308        const int numgrids=1;
    … …  
    387325/*}}}*/
    388326/*FUNCTION Sing::CreatePVector {{{1*/
    389 void  Sing::CreatePVector(Vec pg,void* inputs,int analysis_type,int sub_analysis_type){
     327void  Sing::CreatePVector(Vec pg,int analysis_type,int sub_analysis_type){
    390328       
    391329        /*Just branch to the correct load generator, according to the type of analysis we are carrying out: */
    … …  
    403341/*FUNCTION Sing::CreatePVectorDiagnosticHutter {{{1*/
    404342
    405 void Sing::CreatePVectorDiagnosticHutter( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type){
     343void Sing::CreatePVectorDiagnosticHutter( Vec pg, int analysis_type,int sub_analysis_type){
    406344       
    407345       
    … …  
    561499/*}}}*/
    562500/*FUNCTION Sing::GetOnBed {{{1*/
    563 int   Sing::GetOnBed(){
     501bool   Sing::GetOnBed(){
    564502        ISSMERROR(" not supported yet!");
    565503}

  • issm/trunk/src/c/objects/Sing.h

    r3599 r3612  
    1010#include "./Matice.h"
    1111#include "./Matpar.h"
    12 #include "./ParameterInputs.h"
    13 #include "./ElementProperties.h"
    1412#include "../ModelProcessorx/IoModel.h"
    1513#include "./Hook.h"
    1614
    1715class Hook;
    18 class ElementProperties;
    1916
    2017class Sing: public Element{
    … …  
    3027                Hook hnumpar; //hook to 1 numpar
    3128
    32                 ElementProperties properties;
     29                Inputs* inputs;
    3330
    3431        public:
    … …  
    3633                /*constructors, destructors: {{{1*/
    3734                Sing();
    38                 Sing(int sing_id,int* sing_node_ids, int sing_matice_id, int sing_matpar_id, int sing_numpar_id, ElementProperties* sing_properties);
    39                 Sing(int sing_id,Hook* sing_hnodes, Hook* sing_hmatice, Hook* sing_hmatpar, Hook* sing_hnumpar, ElementProperties* sing_properties);
     35                Sing(int sing_id,int* sing_node_ids, int sing_matice_id, int sing_matpar_id, int sing_numpar_id);
     36                Sing(int sing_id,Hook* sing_hnodes, Hook* sing_hmatice, Hook* sing_hmatpar, Hook* sing_hnumpar, Inputs* sing_inputs);
    4037                Sing(int i, IoModel* iomodel);
    4138                ~Sing();
    … …  
    5653                /*}}}*/
    5754                /*numerics: {{{1*/
    58                 void  CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    59                 void  CreatePVector(Vec pg, void* inputs, int analysis_type,int sub_analysis_type);
    60                 void  UpdateFromInputs(void* inputs);
     55                void  CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type);
     56                void  CreatePVector(Vec pg,  int analysis_type,int sub_analysis_type);
    6157                void  UpdateInputs(double* solution, int analysis_type, int sub_analysis_type);
    6258                void  GetDofList(int* doflist,int* pnumberofdofs);
    6359                void  GetDofList1(int* doflist);
    64                 void  CreateKMatrixDiagnosticHutter(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
     60                void  CreateKMatrixDiagnosticHutter(Mat Kgg,int analysis_type,int sub_analysis_type);
    6561                void  GetParameterValue(double* pp, double* plist, double* gauss_l1l2l3);
    66                 void  CreatePVectorDiagnosticHutter(Vec pg,void* inputs,int analysis_type,int sub_analysis_type);
     62                void  CreatePVectorDiagnosticHutter(Vec pg,int analysis_type,int sub_analysis_type);
    6763                void* GetMatPar();
    68                 void  ComputeBasalStress(Vec sigma_b,void* inputs,int analysis_type,int sub_analysis_type);
    69                 void  ComputePressure(Vec p_g,void* inputs,int analysis_type,int sub_analysis_type);
    70                 void  ComputeStrainRate(Vec eps,void* inputs,int analysis_type,int sub_analysis_type);
     64                void  ComputeBasalStress(Vec sigma_b,int analysis_type,int sub_analysis_type);
     65                void  ComputePressure(Vec p_g,int analysis_type,int sub_analysis_type);
     66                void  ComputeStrainRate(Vec eps,int analysis_type,int sub_analysis_type);
    7167                void  GetNodes(void** vpnodes);
    7268                /*}}}*/
    7369                /*not implemented: {{{1*/
    74                 int   GetShelf();
    75                 int   GetOnBed();
     70                bool   GetShelf();
     71                bool   GetOnBed();
    7672                void  GetBedList(double*);
    7773                void  GetThicknessList(double* thickness_list);

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

    r3599 r3612  
    1919#include "../shared/shared.h"
    2020#include "../DataSet/DataSet.h"
     21#include "../DataSet/Inputs.h"
    2122#include "../include/typedefs.h"
    2223#include "../include/macros.h"
    … …  
    2930}
    3031/*}}}*/
    31 /*FUNCTION Tria::Tria(int id, int* node_ids, int matice_id, int matpar_id, int numpar_id, ElementProperties* properties){{{1*/
    32 Tria::Tria(int tria_id,int* tria_node_ids, int tria_matice_id, int tria_matpar_id, int tria_numpar_id, ElementProperties* triaproperties):
     32/*FUNCTION Tria::Tria(int id, int* node_ids, int matice_id, int matpar_id, int numpar_id){{{1*/
     33Tria::Tria(int tria_id,int* tria_node_ids, int tria_matice_id, int tria_matpar_id, int tria_numpar_id):
    3334        hnodes(tria_node_ids,3),
    3435        hmatice(&tria_matice_id,1),
    3536        hmatpar(&tria_matpar_id,1),
    36         hnumpar(&tria_numpar_id,1),
    37         properties(triaproperties)
     37        hnumpar(&tria_numpar_id,1)
    3838{
    3939
    4040        /*all the initialization has been done by the initializer, just fill in the id: */
    4141        this->id=tria_id;
    42         this->inputs=new DataSet();
     42        this->inputs=new Inputs();
    4343
    4444        return;
    4545}
    4646/*}}}*/
    47 /*FUNCTION Tria::Tria(int id, Hook* hnodes, Hook* hmatice, Hook* hmatpar, Hook* hnumpar, ElementProperties* properties,DataSet* tria_inputs) {{{1*/
    48 Tria::Tria(int tria_id,Hook* tria_hnodes, Hook* tria_hmatice, Hook* tria_hmatpar, Hook* tria_hnumpar, ElementProperties* tria_properties,DataSet* tria_inputs):
     47/*FUNCTION Tria::Tria(int id, Hook* hnodes, Hook* hmatice, Hook* hmatpar, Hook* hnumpar, Inputs* tria_inputs) {{{1*/
     48Tria::Tria(int tria_id,Hook* tria_hnodes, Hook* tria_hmatice, Hook* tria_hmatpar, Hook* tria_hnumpar, Inputs* tria_inputs):
    4949        hnodes(tria_hnodes),
    5050        hmatice(tria_hmatice),
    5151        hmatpar(tria_hmatpar),
    52         hnumpar(tria_hnumpar),
    53         properties(tria_properties)
     52        hnumpar(tria_hnumpar)
    5453{
    5554
    … …  
    5756        this->id=tria_id;
    5857        if(tria_inputs){
    59                 this->inputs=tria_inputs->Copy();
     58                this->inputs=(Inputs*)tria_inputs->Copy();
    6059        }
    6160        else{
    62                 this->inputs=new DataSet();
     61                this->inputs=new Inputs();
    6362        }
    6463
    … …  
    6766/*}}}*/
    6867/*FUNCTION Tria::Tria(int i, IoModel* iomodel){{{1*/
    69 Tria::Tria(int i, IoModel* iomodel){ //i is the element index
    70 
    71         int j;
     68Tria::Tria(int index, IoModel* iomodel){ //i is the element index
     69
     70        int i,j;
    7271        int tria_node_ids[3];
    7372        int tria_matice_id;
    7473        int tria_matpar_id;
    7574        int tria_numpar_id;
     75        double nodeinputs[3];
    7676
    7777        /*id: */
    78         this->id=i+1;
    79         this->inputs=new DataSet();
     78        this->id=index+1;
    8079       
    8180        /*hooks: */
    … …  
    8382        if (iomodel->analysis_type==Prognostic2AnalysisEnum || iomodel->analysis_type==Balancedthickness2AnalysisEnum){
    8483                /*Discontinuous Galerkin*/
    85                 tria_node_ids[0]=3*i+1;
    86                 tria_node_ids[1]=3*i+2;
    87                 tria_node_ids[2]=3*i+3;
     84                tria_node_ids[0]=3*index+1;
     85                tria_node_ids[1]=3*index+2;
     86                tria_node_ids[2]=3*index+3;
    8887        }
    8988        else{
    9089                /*Continuous Galerkin*/
    91                 for(j=0;j<3;j++){
    92                         tria_node_ids[j]=(int)*(iomodel->elements+3*i+j); //ids for vertices are in the elements array from Matlab
     90                for(i=0;i<3;i++){
     91                        tria_node_ids[i]=(int)*(iomodel->elements+3*index+i); //ids for vertices are in the elements array from Matlab
    9392                }
    9493        }
    95         tria_matice_id=i+1; //refers to the corresponding ice material object
     94        tria_matice_id=index+1; //refers to the corresponding ice material object
    9695        tria_matpar_id=iomodel->numberofelements+1; //refers to the constant material parameters object
    9796        tria_numpar_id=1; //refers to numerical parameters object
    … …  
    101100        this->hmatpar.Init(&tria_matpar_id,1);
    102101        this->hnumpar.Init(&tria_numpar_id,1);
    103         this->properties.Init(3,tria_node_ids,this->id,iomodel);
     102
     103        //intialize inputs, and add as many inputs per element as requested:
     104        this->inputs=new Inputs();
     105       
     106        if (iomodel->thickness) {
     107                for(i=0;i<3;i++)nodeinputs[i]=iomodel->thickness[tria_node_ids[i]-1];
     108                this->inputs->AddInput(new TriaVertexInput(ThicknessEnum,nodeinputs));
     109        }
     110        if (iomodel->surface) {
     111                for(i=0;i<3;i++)nodeinputs[i]=iomodel->surface[tria_node_ids[i]-1];
     112                this->inputs->AddInput(new TriaVertexInput(SurfaceEnum,nodeinputs));
     113        }
     114        if (iomodel->bed) {
     115                for(i=0;i<3;i++)nodeinputs[i]=iomodel->bed[tria_node_ids[i]-1];
     116                this->inputs->AddInput(new TriaVertexInput(BedEnum,nodeinputs));
     117        }
     118        if (iomodel->drag_coefficient) {
     119                for(i=0;i<3;i++)nodeinputs[i]=iomodel->drag_coefficient[tria_node_ids[i]-1];
     120                this->inputs->AddInput(new TriaVertexInput(DragCoefficientEnum,nodeinputs));
     121
     122                if (iomodel->drag_p) this->inputs->AddInput(new DoubleInput(DragPEnum,iomodel->drag_p[index]));
     123                if (iomodel->drag_q) this->inputs->AddInput(new DoubleInput(DragQEnum,iomodel->drag_q[index]));
     124                this->inputs->AddInput(new IntInput(DragTypeEnum,iomodel->drag_type));
     125
     126        }
     127        if (iomodel->melting_rate) {
     128                for(i=0;i<3;i++)nodeinputs[i]=iomodel->melting_rate[tria_node_ids[i]-1];
     129                this->inputs->AddInput(new TriaVertexInput(MeltingRateEnum,nodeinputs));
     130        }
     131        if (iomodel->accumulation_rate) {
     132                for(i=0;i<3;i++)nodeinputs[i]=iomodel->accumulation_rate[tria_node_ids[i]-1];
     133                this->inputs->AddInput(new TriaVertexInput(AccumulationRateEnum,nodeinputs));
     134        }
     135        if (iomodel->geothermalflux) {
     136                for(i=0;i<3;i++)nodeinputs[i]=iomodel->geothermalflux[tria_node_ids[i]-1];
     137                this->inputs->AddInput(new TriaVertexInput(GeothermalFluxEnum,nodeinputs));
     138        }       
     139
     140        if (iomodel->elementoniceshelf) this->inputs->AddInput(new BoolInput(ElementOnIceShelfEnum,(IssmBool)iomodel->elementoniceshelf[index]));
     141        if (iomodel->elementonbed) this->inputs->AddInput(new BoolInput(ElementOnBedEnum,(IssmBool)iomodel->elementonbed[index]));
     142        if (iomodel->elementonwater) this->inputs->AddInput(new BoolInput(ElementOnWaterEnum,(IssmBool)iomodel->elementonwater[index]));
     143        if (iomodel->elementonsurface) this->inputs->AddInput(new BoolInput(ElementOnSurfaceEnum,(IssmBool)iomodel->elementonsurface[index]));
     144
    104145
    105146}
    … …  
    141182Object* Tria::copy() {
    142183
    143         return new Tria(this->id,&this->hnodes,&this->hmatice,&this->hmatpar,&this->hnumpar,&this->properties,this->inputs);
     184        return new Tria(this->id,&this->hnodes,&this->hmatice,&this->hmatpar,&this->hnumpar,this->inputs);
    144185
    145186}
    … …  
    166207        hnumpar.Demarshall(&marshalled_dataset);
    167208       
    168         /*demarshall properties: */
    169         properties.Demarshall(&marshalled_dataset);
    170        
    171         inputs=DataSetDemarshallRaw(&marshalled_dataset);
     209        /*demarshall inputs: */
     210        inputs=(Inputs*)DataSetDemarshallRaw(&marshalled_dataset);
    172211
    173212        /*return: */
    … …  
    186225        hmatpar.DeepEcho();
    187226        hnumpar.DeepEcho();
    188         properties.DeepEcho();
    189         printf("Element inputs\n");
     227        printf("   inputs\n");
    190228        inputs->DeepEcho();
    191229
    … …  
    203241        hmatpar.Echo();
    204242        hnumpar.Echo();
    205         properties.Echo();
    206         printf("Element inputs\n");
     243        printf("   inputs\n");
    207244        inputs->Echo();
    208245
    … …  
    236273        hnumpar.Marshall(&marshalled_dataset);
    237274
    238         /*Marshall properties: */
    239         properties.Marshall(&marshalled_dataset);
    240 
    241275        /*Marshall inputs: */
    242276        marshalled_inputs_size=inputs->MarshallSize();
    … …  
    257291                +hmatpar.MarshallSize()
    258292                +hnumpar.MarshallSize()
    259                 +properties.MarshallSize()
    260293                +inputs->MarshallSize()
    261294                +sizeof(int); //sizeof(int) for enum type
    … …  
    264297
    265298/*Updates: */
    266 /*FUNCTION Tria::UpdateFromInputs {{{1*/
    267 void  Tria::UpdateFromInputs(void* vinputs){
     299/*FUNCTION Tria::UpdateFromDakota {{{1*/
     300void  Tria::UpdateFromDakota(void* vinputs){
    268301
    269302        int     i;
    … …  
    280313        Matice* matice=NULL;
    281314        Numpar* numpar=NULL;
    282 
    283         ParameterInputs* inputs=NULL;
    284 
    285         /*recover pointers: */
    286         inputs=(ParameterInputs*)vinputs;
    287 
    288         /*recover objects from hooks: */
    289         nodes=(Node**)hnodes.deliverp();
    290         matpar=(Matpar*)hmatpar.delivers();
    291         matice=(Matice*)hmatice.delivers();
    292         numpar=(Numpar*)hnumpar.delivers();
    293 
    294         /*Update internal data if inputs holds new values: */
    295         //if (id==1) printf("WARNING if QMU: no hydrostatic equilibrium is applied here (conflict with prognostic, which does not have matpar)\n");
    296         //For now
    297         if(this->properties.h) inputs->Recover("thickness",&this->properties.h[0],1,dofs,3,(void**)nodes);
    298         //Later
    299         /*
    300                 if(inputs->Recover("thickness",&new_h[0],1,dofs,3,(void**)nodes)){
    301         //density, needed later:
    302         double di=(this->matpar->GetRhoIce()/this->matpar->GetRhoWater());
    303         //Go through grids:
    304         for (i=0;i<3;i++){
    305         if(nodes[i]->IsOnShelf()){
    306         this->b[i]=this->b[i]-di*(new_h[i]-h[i]); //hydrostatic equilibrium;
    307         }
    308         this->s[i]=this->b[i]+new_h[i];
    309         this->h[i]=new_h[i];
    310         }
    311         }
    312         */
    313         if (this->properties.s)              inputs->Recover("surface",&this->properties.s[0],1,dofs,3,(void**)nodes);
    314         if (this->properties.b)              inputs->Recover("bed",&this->properties.b[0],1,dofs,3,(void**)nodes);
    315         if (this->properties.k)              inputs->Recover("drag",&this->properties.k[0],1,dofs,3,(void**)nodes);
    316         if (this->properties.melting)        inputs->Recover("melting",&this->properties.melting[0],1,dofs,3,(void**)nodes);
    317         if (this->properties.accumulation)   inputs->Recover("accumulation",&this->properties.accumulation[0],1,dofs,3,(void**)nodes);
    318         if (this->properties.geothermalflux) inputs->Recover("geothermalflux",&this->properties.geothermalflux[0],1,dofs,3,(void**)nodes);
    319        
    320         //Update material if necessary
    321         if(inputs->Recover("temperature_average",&temperature_list[0],1,dofs,3,(void**)nodes)){
    322                 temperature_average=(temperature_list[0]+temperature_list[1]+temperature_list[2])/3.0;
    323                 B_average=Paterson(temperature_average);
    324                 matice->SetB(B_average);
    325         }
    326        
    327         if(inputs->Recover("B",&B_list[0],1,dofs,3,(void**)nodes)){
    328                 B_average=(B_list[0]+B_list[1]+B_list[2])/3.0;
    329                 matice->SetB(B_average);
    330         }
    331 
    332 }
    333 /*}}}*/
    334 /*FUNCTION Tria::UpdateFromDakota {{{1*/
    335 void  Tria::UpdateFromDakota(void* vinputs){
    336 
    337         int     i;
    338         int     dofs[1]={0};
    339         double  temperature_list[3];
    340         double  temperature_average;
    341         double  B_list[3];
    342         double  B_average;
    343         double  new_h[3];
    344 
    345         /*dynamic objects pointed to by hooks: */
    346         Node**  nodes=NULL;
    347         Matpar* matpar=NULL;
    348         Matice* matice=NULL;
    349         Numpar* numpar=NULL;
    350 
    351         ParameterInputs* inputs=NULL;
    352 
    353         /*recover pointers: */
    354         inputs=(ParameterInputs*)vinputs;
    355315
    356316        /*recover objects from hooks: */
    … …  
    458418
    459419        /*Add vx and vy as inputs to the tria element: */
    460         this->inputs->AddEinput(new TriaVertexInput(VxEnum,vx));
    461         this->inputs->AddEinput(new TriaVertexInput(VyEnum,vy));
     420        this->inputs->AddInput(new TriaVertexInput(VxEnum,vx));
     421        this->inputs->AddInput(new TriaVertexInput(VyEnum,vy));
    462422}
    463423
    … …  
    497457/*Object functions*/
    498458/*FUNCTION Tria::ComputeBasalStress {{{1*/
    499 void  Tria::ComputeBasalStress(Vec eps,void* vinputs,int analysis_type,int sub_analysis_type){
     459void  Tria::ComputeBasalStress(Vec eps,int analysis_type,int sub_analysis_type){
    500460
    501461        int i;
    … …  
    523483/*}}}*/
    524484/*FUNCTION Tria::ComputePressure {{{1*/
    525 void  Tria::ComputePressure(Vec pg,void* vinputs,int analysis_type,int sub_analysis_type){
     485void  Tria::ComputePressure(Vec pg,int analysis_type,int sub_analysis_type){
    526486
    527487        int i;
    528         const int numgrids=3;
    529         int doflist[numgrids];
    530         double pressure[numgrids];
     488        const int numvertices=3;
     489        int doflist[numvertices];
     490        double pressure[numvertices];
     491        double thickness[numvertices];
    531492        double rho_ice,g;
     493        double gauss[numvertices][numvertices]={{1,0,0},{0,1,0},{0,0,1}};
    532494       
    533495        /*dynamic objects pointed to by hooks: */
    … …  
    537499        Numpar* numpar=NULL;
    538500
    539        
    540         /*Get dof list on which we will plug the pressure values: */
    541         GetDofList1(&doflist[0]);
    542 
    543501        /*recover objects from hooks: */
    544502        nodes=(Node**)hnodes.deliverp();
    … …  
    547505        numpar=(Numpar*)hnumpar.delivers();
    548506
     507        /*Get dof list on which we will plug the pressure values: */
     508        GetDofList1(&doflist[0]);
     509
    549510        /*pressure is lithostatic: */
    550511        rho_ice=matpar->GetRhoIce();
    551512        g=matpar->GetG();
    552         for(i=0;i<numgrids;i++){
    553                 pressure[i]=rho_ice*g*this->properties.h[i];
     513
     514        /*recover value of thickness at gauss points (0,0,1), (0,1,0),(1,0,0): */
     515        inputs->GetParameterValues(&thickness[0],&gauss[0][0],3,ThicknessEnum);
     516
     517        for(i=0;i<numvertices;i++){
     518                pressure[i]=rho_ice*g*thickness[i];
    554519        }
    555520
    556521        /*plug local pressure values into global pressure vector: */
    557         VecSetValues(pg,numgrids,doflist,(const double*)pressure,INSERT_VALUES);
     522        VecSetValues(pg,numvertices,doflist,(const double*)pressure,INSERT_VALUES);
    558523
    559524}
    560525/*}}}*/
    561526/*FUNCTION Tria::ComputeStrainRate {{{1*/
    562 void  Tria::ComputeStrainRate(Vec eps,void* vinputs,int analysis_type,int sub_analysis_type){
     527void  Tria::ComputeStrainRate(Vec eps,int analysis_type,int sub_analysis_type){
    563528
    564529        int i;
    … …  
    586551/*}}}*/
    587552/*FUNCTION Tria::CostFunction {{{1*/
    588 double Tria::CostFunction(void* vinputs,int analysis_type,int sub_analysis_type){
     553double Tria::CostFunction(int analysis_type,int sub_analysis_type){
    589554
    590555        int i;
    … …  
    628593        Numpar* numpar=NULL;
    629594
    630         ParameterInputs* inputs=NULL;
     595        /*inputs: */
     596        bool onwater;
     597        bool shelf;
     598
     599        /*retrieve inputs :*/
     600        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     601        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
    631602
    632603        /*If on water, return 0: */
    633         if(this->properties.onwater)return 0;
    634 
    635         /*recover pointers: */
    636         inputs=(ParameterInputs*)vinputs;
     604        if(onwater)return 0;
    637605
    638606        /*recover objects from hooks: */
    … …  
    646614
    647615        /*First, get Misfit*/
    648         Jelem=Misfit(inputs,analysis_type,sub_analysis_type);
     616        Jelem=Misfit(analysis_type,sub_analysis_type);
    649617
    650618          /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
    … …  
    664632                /*Add Tikhonov regularization term to misfit*/
    665633                if (strcmp(numpar->control_type,"drag")==0){
    666                         if (!this->properties.shelf){
    667 
    668                                 GetParameterDerivativeValue(&dk[0], &this->properties.k[0],&xyz_list[0][0], gauss_l1l2l3);
     634                        if (shelf){
     635                                inputs->GetParameterDerivativeValue(&dk[0],&xyz_list[0][0],&gauss_l1l2l3[0],DragCoefficientEnum);
    669636                                Jelem+=numpar->cm_noisedmp*1/2*(pow(dk[0],2)+pow(dk[1],2))*Jdet*gauss_weight;
    670637
    … …  
    672639                }
    673640                else if (strcmp(numpar->control_type,"B")==0){
    674 
    675                         if(!inputs->Recover("B",&B[0],1,dofs1,numgrids,(void**)nodes)){
    676                                 ISSMERROR("parameter B not found in input");
    677                         }
    678                         GetParameterDerivativeValue(&dB[0], &B[0],&xyz_list[0][0], gauss_l1l2l3);
     641                        inputs->GetParameterDerivativeValue(&dB[0], &xyz_list[0][0], &gauss_l1l2l3[0],RheologyBEnum);
    679642                        Jelem+=numpar->cm_noisedmp*1/2*(pow(dB[0],2)+pow(dB[1],2))*Jdet*gauss_weight;
    680 
    681643                }
    682644                else{
    … …  
    697659/*FUNCTION Tria::CreateKMatrix {{{1*/
    698660
    699 void  Tria::CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     661void  Tria::CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type){
    700662
    701663        /*Just branch to the correct element stiffness matrix generator, according to the type of analysis we are carrying out: */
    702664        if (analysis_type==ControlAnalysisEnum){
    703665               
    704                 CreateKMatrixDiagnosticHoriz( Kgg,inputs,analysis_type,sub_analysis_type);
     666                CreateKMatrixDiagnosticHoriz( Kgg,analysis_type,sub_analysis_type);
    705667        }
    706668        else if (analysis_type==DiagnosticAnalysisEnum){
    … …  
    708670                if (sub_analysis_type==HorizAnalysisEnum){
    709671
    710                         CreateKMatrixDiagnosticHoriz( Kgg,inputs,analysis_type,sub_analysis_type);
     672                        CreateKMatrixDiagnosticHoriz( Kgg,analysis_type,sub_analysis_type);
    711673                }
    712674                else ISSMERROR("%s%i%s\n","sub_analysis: ",sub_analysis_type," not supported yet");
    … …  
    715677        else if (analysis_type==SlopecomputeAnalysisEnum){
    716678
    717                 CreateKMatrixSlopeCompute( Kgg,inputs,analysis_type,sub_analysis_type);
     679                CreateKMatrixSlopeCompute( Kgg,analysis_type,sub_analysis_type);
    718680        }
    719681        else if (analysis_type==PrognosticAnalysisEnum){
    720682
    721                 CreateKMatrixPrognostic( Kgg,inputs,analysis_type,sub_analysis_type);
     683                CreateKMatrixPrognostic( Kgg,analysis_type,sub_analysis_type);
    722684        }
    723685        else if (analysis_type==Prognostic2AnalysisEnum){
    724686
    725                 CreateKMatrixPrognostic2(Kgg,inputs,analysis_type,sub_analysis_type);
     687                CreateKMatrixPrognostic2(Kgg,analysis_type,sub_analysis_type);
    726688        }
    727689        else if (analysis_type==BalancedthicknessAnalysisEnum){
    728690
    729                 CreateKMatrixBalancedthickness( Kgg,inputs,analysis_type,sub_analysis_type);
     691                CreateKMatrixBalancedthickness( Kgg,analysis_type,sub_analysis_type);
    730692        }
    731693        else if (analysis_type==Balancedthickness2AnalysisEnum){
    732694
    733                 CreateKMatrixBalancedthickness2( Kgg,inputs,analysis_type,sub_analysis_type);
     695                CreateKMatrixBalancedthickness2( Kgg,analysis_type,sub_analysis_type);
    734696        }
    735697        else if (analysis_type==BalancedvelocitiesAnalysisEnum){
    736698
    737                 CreateKMatrixBalancedvelocities( Kgg,inputs,analysis_type,sub_analysis_type);
     699                CreateKMatrixBalancedvelocities( Kgg,analysis_type,sub_analysis_type);
    738700        }
    739701        else{
    … …  
    745707/*}}}*/
    746708/*FUNCTION Tria::CreateKMatrixBalancedthickness {{{1*/
    747 void  Tria::CreateKMatrixBalancedthickness(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     709void  Tria::CreateKMatrixBalancedthickness(Mat Kgg,int analysis_type,int sub_analysis_type){
    748710
    749711        /* local declarations */
    … …  
    766728        double  gauss_weight;
    767729        double  gauss_l1l2l3[3];
     730        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
    768731
    769732        /* matrices: */
    … …  
    782745
    783746        /*input parameters for structural analysis (diagnostic): */
    784         double  vxvy_list[numgrids][2]={0.0};
    785         double  vx_list[numgrids]={0.0};
    786         double  vy_list[numgrids]={0.0};
     747        double  vx_list[numgrids];
     748        double  vy_list[numgrids];
    787749        double  dvx[2];
    788750        double  dvy[2];
    … …  
    790752        double  dvxdx,dvydy;
    791753        double  v_gauss[2]={0.0};
     754
     755
    792756        double  K[2][2]={0.0};
    793757        double  KDL[2][2]={0.0};
    … …  
    801765        Numpar* numpar=NULL;
    802766
    803         ParameterInputs* inputs=NULL;
    804 
    805         /*recover pointers: */
    806         inputs=(ParameterInputs*)vinputs;
    807767
    808768        /*recover objects from hooks: */
    … …  
    812772        numpar=(Numpar*)hnumpar.delivers();
    813773
    814         /*recover extra inputs from users, at current convergence iteration: */
    815         found=inputs->Recover("velocity_average",&vxvy_list[0][0],2,dofs,numgrids,(void**)nodes);
    816         if(!found)ISSMERROR(" could not find velocity_average  in inputs!");
    817 
    818         for(i=0;i<numgrids;i++){
    819                 vx_list[i]=vxvy_list[i][0];
    820                 vy_list[i]=vxvy_list[i][1];
    821         }
    822 
    823774        /* Get node coordinates and dof list: */
    824775        GetVerticesCoordinates(&xyz_list[0][0], nodes, numgrids);
    825776        GetDofList(&doflist[0],&numberofdofspernode);
     777
     778        /*Recover velocity: */
     779        inputs->GetParameterValues(&vx_list[0],&gaussgrids[0][0],3,VxAverageEnum);
     780        inputs->GetParameterValues(&vy_list[0],&gaussgrids[0][0],3,VyAverageEnum);
    826781
    827782        //Create Artificial diffusivity once for all if requested
    … …  
    832787
    833788                //Build K matrix (artificial diffusivity matrix)
    834                 v_gauss[0]=ONETHIRD*(vxvy_list[0][0]+vxvy_list[1][0]+vxvy_list[2][0]);
    835                 v_gauss[1]=ONETHIRD*(vxvy_list[0][1]+vxvy_list[1][1]+vxvy_list[2][1]);
     789                v_gauss[0]=ONETHIRD*(vx_list[0]+vx_list[1]+vx_list[2]);
     790                v_gauss[1]=ONETHIRD*(vy_list[0]+vy_list[1]+vy_list[2]);
    836791
    837792                K[0][0]=pow(Jdettria,(double).5)/2.0*fabs(v_gauss[0]);
    … …  
    858813
    859814                //Get vx, vy and their derivatives at gauss point
    860                 GetParameterValue(&vx, &vx_list[0],gauss_l1l2l3);
    861                 GetParameterValue(&vy, &vy_list[0],gauss_l1l2l3);
    862 
    863                 GetParameterDerivativeValue(&dvx[0], &vx_list[0],&xyz_list[0][0], gauss_l1l2l3);
    864                 GetParameterDerivativeValue(&dvy[0], &vy_list[0],&xyz_list[0][0], gauss_l1l2l3);
     815                this->GetParameterValue(&vx, &vx_list[0],gauss_l1l2l3);
     816                this->GetParameterValue(&vy, &vy_list[0],gauss_l1l2l3);
     817
     818                this->GetParameterDerivativeValue(&dvx[0], &vx_list[0],&xyz_list[0][0], gauss_l1l2l3);
     819                this->GetParameterDerivativeValue(&dvy[0], &vy_list[0],&xyz_list[0][0], gauss_l1l2l3);
    865820
    866821                dvxdx=dvx[0];
    … …  
    922877/*}}}*/
    923878/*FUNCTION Tria::CreateKMatrixBalancedthickness2 {{{1*/
    924 void  Tria::CreateKMatrixBalancedthickness2(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     879void  Tria::CreateKMatrixBalancedthickness2(Mat Kgg,int analysis_type,int sub_analysis_type){
    925880
    926881        /* local declarations */
    … …  
    955910
    956911        /*input parameters for structural analysis (diagnostic): */
    957         double  vx_list[numgrids]={0.0};
    958         double  vy_list[numgrids]={0.0};
    959912        double  vx,vy;
    960913        int     dofs[1]={0};
    … …  
    967920        Numpar* numpar=NULL;
    968921
    969         ParameterInputs* inputs=NULL;
    970 
    971         /*recover pointers: */
    972         inputs=(ParameterInputs*)vinputs;
    973922
    974923        /*recover objects from hooks: */
    … …  
    977926        matice=(Matice*)hmatice.delivers();
    978927        numpar=(Numpar*)hnumpar.delivers();
    979 
    980         /*recover extra inputs from users, at current convergence iteration: */
    981         found=inputs->Recover("vx_average",&vx_list[0],1,dofs,numgrids,(void**)nodes);
    982         if(!found)ISSMERROR(" could not find vx_average in inputs!");
    983         found=inputs->Recover("vy_average",&vy_list[0],1,dofs,numgrids,(void**)nodes);
    984         if(!found)ISSMERROR(" could not find vy_average in inputs!");
    985928
    986929        /* Get node coordinates and dof list: */
    … …  
    1008951
    1009952                //Get vx, vy and their derivatives at gauss point
    1010                 GetParameterValue(&vx, &vx_list[0],gauss_l1l2l3);
    1011                 GetParameterValue(&vy, &vy_list[0],gauss_l1l2l3);
     953                inputs->GetParameterValue(&vx, &gauss_l1l2l3[0],VxEnum);
     954                inputs->GetParameterValue(&vy, &gauss_l1l2l3[0],VyEnum);
    1012955
    1013956                DL_scalar=-gauss_weight*Jdettria;
    … …  
    1039982/*}}}*/
    1040983/*FUNCTION Tria::CreateKMatrixBalancedvelocities {{{1*/
    1041 void  Tria::CreateKMatrixBalancedvelocities(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     984void  Tria::CreateKMatrixBalancedvelocities(Mat Kgg,int analysis_type,int sub_analysis_type){
    1042985
    1043986        /* local declarations */
    … …  
    10601003        double  gauss_weight;
    10611004        double  gauss_l1l2l3[3];
     1005        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
    10621006
    10631007        /* matrices: */
    … …  
    10761020        /*input parameters for structural analysis (diagnostic): */
    10771021        double  surface_normal[3];
     1022        double  surface_list[3];
    10781023        double  nx,ny,norm;
    1079         double  vxvy_list[numgrids][2]={0.0};
    10801024        double  vx_list[numgrids]={0.0};
    10811025        double  vy_list[numgrids]={0.0};
    … …  
    10961040        Numpar* numpar=NULL;
    10971041
    1098         ParameterInputs* inputs=NULL;
    1099 
    1100         /*recover pointers: */
    1101         inputs=(ParameterInputs*)vinputs;
    1102 
    11031042        /*recover objects from hooks: */
    11041043        nodes=(Node**)hnodes.deliverp();
    … …  
    11071046        numpar=(Numpar*)hnumpar.delivers();
    11081047
    1109         /*recover extra inputs from users, at current convergence iteration: */
    1110         found=inputs->Recover("velocity_average",&vxvy_list[0][0],2,dofs,numgrids,(void**)nodes);
    1111         if(!found)ISSMERROR(" could not find velocity_average  in inputs!");
    1112 
    1113         for(i=0;i<numgrids;i++){
    1114                 vx_list[i]=vxvy_list[i][0];
    1115                 vy_list[i]=vxvy_list[i][1];
    1116         }
     1048        /*Recover velocity: */
     1049        inputs->GetParameterValues(&vx_list[0],&gaussgrids[0][0],3,VxAverageEnum);
     1050        inputs->GetParameterValues(&vy_list[0],&gaussgrids[0][0],3,VyAverageEnum);
    11171051
    11181052        /* Get node coordinates and dof list: */
    … …  
    11211055
    11221056        /*Modify z so that it reflects the surface*/
    1123         for(i=0;i<numgrids;i++) xyz_list[i][2]=this->properties.s[i];
     1057        inputs->GetParameterValues(&surface_list[0],&gaussgrids[0][0],3,SurfaceEnum);
     1058        for(i=0;i<numgrids;i++) xyz_list[i][2]=surface_list[i];
    11241059
    11251060        /*Get normal vector to the surface*/
    … …  
    11461081
    11471082                //Build K matrix (artificial diffusivity matrix)
    1148                 v_gauss[0]=ONETHIRD*(vxvy_list[0][0]+vxvy_list[1][0]+vxvy_list[2][0]);
    1149                 v_gauss[1]=ONETHIRD*(vxvy_list[0][1]+vxvy_list[1][1]+vxvy_list[2][1]);
     1083                v_gauss[0]=ONETHIRD*(vx_list[0]+vx_list[1]+vx_list[2]);
     1084                v_gauss[1]=ONETHIRD*(vy_list[0]+vy_list[1]+vy_list[2]);
    11501085
    11511086                K[0][0]=pow(10,2)*pow(Jdettria,(double).5)/2.0*fabs(v_gauss[0]); //pow should be zero!!
    … …  
    12271162/*}}}*/
    12281163/*FUNCTION Tria::CreateKMatrixDiagnosticHoriz {{{1*/
    1229 void  Tria::CreateKMatrixDiagnosticHoriz(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     1164void  Tria::CreateKMatrixDiagnosticHoriz(Mat Kgg,int analysis_type,int sub_analysis_type){
    12301165
    12311166        /* local declarations */
    … …  
    12701205
    12711206        /*input parameters for structural analysis (diagnostic): */
    1272         double  vxvy_list[numgrids][2]={{0,0},{0,0},{0,0}};
    1273         double  oldvxvy_list[numgrids][2]={{0,0},{0,0},{0,0}};
    12741207        double  thickness;
    12751208        int     dofs[2]={0,1};
    … …  
    12811214        Numpar* numpar=NULL;
    12821215
    1283         ParameterInputs* inputs=NULL;
     1216        /*inputs: */
     1217        bool onwater,shelf;
     1218
     1219        /*retrieve inputs :*/
     1220        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     1221        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
    12841222
    12851223        /*First, if we are on water, return empty matrix: */
    1286         if(this->properties.onwater) return;
    1287 
    1288         /*recover pointers: */
    1289         inputs=(ParameterInputs*)vinputs;
     1224        if(onwater) return;
    12901225
    12911226        /*recover objects from hooks: */
    … …  
    12941229        matice=(Matice*)hmatice.delivers();
    12951230        numpar=(Numpar*)hnumpar.delivers();
    1296 
    1297         /*recover extra inputs from users, at current convergence iteration: */
    1298         inputs->Recover("velocity",&vxvy_list[0][0],2,dofs,numgrids,(void**)nodes);
    1299         inputs->Recover("old_velocity",&oldvxvy_list[0][0],2,dofs,numgrids,(void**)nodes);
    13001231
    13011232        /* Get node coordinates and dof list: */
    … …  
    13191250
    13201251                /*Compute thickness at gaussian point: */
    1321                 GetParameterValue(&thickness, &this->properties.h[0],gauss_l1l2l3);
     1252                inputs->GetParameterValue(&thickness, gauss_l1l2l3,ThicknessEnum);
    13221253
    13231254                /*Get strain rate from velocity: */
    1324                 GetStrainRate(&epsilon[0],&vxvy_list[0][0],&xyz_list[0][0],gauss_l1l2l3);
    1325                 GetStrainRate(&oldepsilon[0],&oldvxvy_list[0][0],&xyz_list[0][0],gauss_l1l2l3);
     1255                inputs->GetStrainRate(&epsilon[0],&xyz_list[0][0],gauss_l1l2l3,VxEnum,VyEnum);
     1256                inputs->GetStrainRate(&oldepsilon[0],&xyz_list[0][0],gauss_l1l2l3,VxOldEnum,VyOldEnum);
    13261257
    13271258                /*Get viscosity: */
    … …  
    13601291
    13611292        /*Do not forget to include friction: */
    1362         if(!this->properties.shelf){
    1363                 CreateKMatrixDiagnosticHorizFriction(Kgg,inputs,analysis_type,sub_analysis_type);
     1293        if(!shelf){
     1294                CreateKMatrixDiagnosticHorizFriction(Kgg,analysis_type,sub_analysis_type);
    13641295        }
    13651296
    … …  
    13731304/*}}}*/
    13741305/*FUNCTION Tria::CreateKMatrixDiagnosticHorizFriction {{{1*/
    1375 void  Tria::CreateKMatrixDiagnosticHorizFriction(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     1306void  Tria::CreateKMatrixDiagnosticHorizFriction(Mat Kgg,int analysis_type,int sub_analysis_type){
    13761307
    13771308
    … …  
    13941325        double  gauss_weight;
    13951326        double  gauss_l1l2l3[3];
     1327        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
    13961328
    13971329        /* matrices: */
    … …  
    14111343
    14121344        /*input parameters for structural analysis (diagnostic): */
    1413         double  vxvy_list[numgrids][2]={{0,0},{0,0},{0,0}};
    1414         int     dofs[2]={0,1};
     1345        double  vx_list[numgrids];
     1346        double  vy_list[numgrids];
     1347        double  thickness_list[numgrids];
     1348        double  bed_list[numgrids];
     1349        double  dragcoefficient_list[numgrids];
     1350        double  drag_p,drag_q;
    14151351
    14161352        /*friction: */
    … …  
    14201356        double MAXSLOPE=.06; // 6 %
    14211357        double MOUNTAINKEXPONENT=10;
     1358
     1359        /*inputs: */
     1360        bool shelf;
     1361        int  drag_type;
    14221362
    14231363        /*dynamic objects pointed to by hooks: */
    … …  
    14271367        Numpar* numpar=NULL;
    14281368
    1429         ParameterInputs* inputs=NULL;
    1430 
    1431         /*recover pointers: */
    1432         inputs=(ParameterInputs*)vinputs;
    14331369
    14341370        /*recover objects from hooks: */
    … …  
    14371373        matice=(Matice*)hmatice.delivers();
    14381374        numpar=(Numpar*)hnumpar.delivers();
     1375
     1376        /*retrieve inputs :*/
     1377        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
     1378        inputs->GetParameterValue(&drag_type,DragTypeEnum);
    14391379       
    14401380        /* Get node coordinates and dof list: */
    … …  
    14451385        for(i=0;i<numdof;i++) for(j=0;j<numdof;j++) Ke_gg[i][j]=0.0;
    14461386
    1447         if (this->properties.shelf){
     1387        if (shelf){
    14481388                /*no friction, do nothing*/
    14491389                return;
    14501390        }
    14511391
    1452         if (this->properties.friction_type!=2)ISSMERROR(" non-viscous friction not supported yet!");
    1453 
    1454         /*recover extra inputs from users, at current convergence iteration: */
    1455         inputs->Recover("velocity",&vxvy_list[0][0],2,dofs,numgrids,(void**)nodes);
     1392        if (drag_type!=2)ISSMERROR(" non-viscous friction not supported yet!");
     1393
     1394        /*Recover inputs: */
     1395        inputs->GetParameterValues(&vx_list[0],&gaussgrids[0][0],3,VxAverageEnum);
     1396        inputs->GetParameterValues(&vy_list[0],&gaussgrids[0][0],3,VyAverageEnum);
     1397        inputs->GetParameterValues(&dragcoefficient_list[0],&gaussgrids[0][0],3,DragCoefficientEnum);
     1398        inputs->GetParameterValues(&bed_list[0],&gaussgrids[0][0],3,BedEnum);
     1399        inputs->GetParameterValues(&thickness_list[0],&gaussgrids[0][0],3,ThicknessEnum);
     1400        inputs->GetParameterValue(&drag_p,DragPEnum);
     1401        inputs->GetParameterValue(&drag_q,DragQEnum);
    14561402
    14571403        /*Build alpha2_list used by drag stiffness matrix*/
    … …  
    14651411        friction->rho_ice=matpar->GetRhoIce();
    14661412        friction->rho_water=matpar->GetRhoWater();
    1467         friction->K=&this->properties.k[0];
    1468         friction->bed=&this->properties.b[0];
    1469         friction->thickness=&this->properties.h[0];
    1470         friction->velocities=&vxvy_list[0][0];
    1471         friction->p=this->properties.p;
    1472         friction->q=this->properties.q;
     1413        friction->K=&dragcoefficient_list[0];
     1414        friction->bed=&bed_list[0];
     1415        friction->thickness=&thickness_list[0];
     1416        friction->vx=&vx_list[0];
     1417        friction->vy=&vy_list[0];
     1418        friction->p=drag_p;
     1419        friction->q=drag_q;
    14731420
    14741421        /*Compute alpha2_list: */
    … …  
    14921439                // If we have a slope > 6% for this element,  it means  we are on a mountain. In this particular case,
    14931440                //velocity should be = 0. To achieve this result, we set alpha2_list to a very high value: */
    1494                 GetParameterDerivativeValue(&slope[0], &this->properties.s[0],&xyz_list[0][0], gauss_l1l2l3);
     1441                inputs->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],&gauss_l1l2l3[0],SurfaceEnum);
    14951442                slope_magnitude=sqrt(pow(slope[0],2)+pow(slope[1],2));
    14961443
    … …  
    15371484/*}}}*/
    15381485/*FUNCTION Tria::CreateKMatrixDiagnosticSurfaceVert {{{1*/
    1539 void  Tria::CreateKMatrixDiagnosticSurfaceVert(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     1486void  Tria::CreateKMatrixDiagnosticSurfaceVert(Mat Kgg,int analysis_type,int sub_analysis_type){
    15401487
    15411488        int i,j;
    … …  
    15841531        Numpar* numpar=NULL;
    15851532
    1586         ParameterInputs* inputs=NULL;
    1587 
    1588         /*recover pointers: */
    1589         inputs=(ParameterInputs*)vinputs;
    1590 
    15911533        /*recover objects from hooks: */
    15921534        nodes=(Node**)hnodes.deliverp();
    … …  
    16751617/*}}}*/
    16761618/*FUNCTION Tria::CreateKMatrixMelting {{{1*/
    1677 void  Tria::CreateKMatrixMelting(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     1619void  Tria::CreateKMatrixMelting(Mat Kgg,int analysis_type,int sub_analysis_type){
    16781620
    16791621        /*indexing: */
    … …  
    17711713/*}}}*/
    17721714/*FUNCTION Tria::CreateKMatrixPrognostic {{{1*/
    1773 void  Tria::CreateKMatrixPrognostic(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     1715void  Tria::CreateKMatrixPrognostic(Mat Kgg,int analysis_type,int sub_analysis_type){
    17741716
    17751717        /* local declarations */
    … …  
    17921734        double  gauss_weight;
    17931735        double  gauss_l1l2l3[3];
     1736        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
    17941737
    17951738        /* matrices: */
    … …  
    18071750
    18081751        /*input parameters for structural analysis (diagnostic): */
    1809         double  vxvy_list[numgrids][2]={0.0};
    18101752        double  vx_list[numgrids]={0.0};
    18111753        double  vy_list[numgrids]={0.0};
    … …  
    18271769        Numpar* numpar=NULL;
    18281770
    1829         ParameterInputs* inputs=NULL;
    1830 
    1831         /*recover pointers: */
    1832         inputs=(ParameterInputs*)vinputs;
    1833 
    18341771        /*recover objects from hooks: */
    18351772        nodes=(Node**)hnodes.deliverp();
    … …  
    18381775        numpar=(Numpar*)hnumpar.delivers();
    18391776
    1840         /*recover extra inputs from users, at current convergence iteration: */
    1841         found=inputs->Recover("velocity_average",&vxvy_list[0][0],2,dofs,numgrids,(void**)nodes);
    1842         if(!found)ISSMERROR(" could not find velocity_average  in inputs!");
    1843 
    1844         for(i=0;i<numgrids;i++){
    1845                 vx_list[i]=vxvy_list[i][0];
    1846                 vy_list[i]=vxvy_list[i][1];
    1847         }
    1848 
    1849         found=inputs->Recover("dt",&dt);
    1850         if(!found)ISSMERROR(" could not find dt in inputs!");
     1777        /*recover dt: */
     1778        dt=numpar->dt;
    18511779
    18521780        /* Get node coordinates and dof list: */
    18531781        GetVerticesCoordinates(&xyz_list[0][0], nodes, numgrids);
    18541782        GetDofList(&doflist[0],&numberofdofspernode);
     1783
     1784        /*Recover velocity: */
     1785        inputs->GetParameterValues(&vx_list[0],&gaussgrids[0][0],3,VxAverageEnum);
     1786        inputs->GetParameterValues(&vy_list[0],&gaussgrids[0][0],3,VyAverageEnum);
    18551787
    18561788        //Create Artificial diffusivity once for all if requested
    … …  
    18611793
    18621794                //Build K matrix (artificial diffusivity matrix)
    1863                 v_gauss[0]=ONETHIRD*(vxvy_list[0][0]+vxvy_list[1][0]+vxvy_list[2][0]);
    1864                 v_gauss[1]=ONETHIRD*(vxvy_list[0][1]+vxvy_list[1][1]+vxvy_list[2][1]);
     1795                v_gauss[0]=ONETHIRD*(vx_list[0]+vx_list[1]+vx_list[2]);
     1796                v_gauss[1]=ONETHIRD*(vy_list[0]+vy_list[1]+vy_list[2]);
    18651797
    18661798                K[0][0]=pow(Jdettria,(double).5)/2.0*fabs(v_gauss[0]);
    … …  
    19661898/*}}}*/
    19671899/*FUNCTION Tria::CreateKMatrixPrognostic2 {{{1*/
    1968 void  Tria::CreateKMatrixPrognostic2(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     1900void  Tria::CreateKMatrixPrognostic2(Mat Kgg,int analysis_type,int sub_analysis_type){
    19691901
    19701902        /* local declarations */
    … …  
    20011933
    20021934        /*input parameters for structural analysis (diagnostic): */
    2003         double  vx_list[numgrids]={0.0};
    2004         double  vy_list[numgrids]={0.0};
    20051935        double  vx,vy;
    20061936        double  dt;
    … …  
    20141944        Numpar* numpar=NULL;
    20151945
    2016         ParameterInputs* inputs=NULL;
    2017 
    2018         /*recover pointers: */
    2019         inputs=(ParameterInputs*)vinputs;
    2020 
    20211946        /*recover objects from hooks: */
    20221947        nodes=(Node**)   hnodes.deliverp();
    … …  
    20251950        numpar=(Numpar*)hnumpar.delivers();
    20261951
    2027         /*recover extra inputs from users, at current convergence iteration: */
    2028         found=inputs->Recover("vx_average",&vx_list[0],1,dofs,numgrids,(void**)nodes);
    2029         if(!found)ISSMERROR(" could not find vx_average in inputs!");
    2030         found=inputs->Recover("vy_average",&vy_list[0],1,dofs,numgrids,(void**)nodes);
    2031         if(!found)ISSMERROR(" could not find vy_average in inputs!");
    2032         found=inputs->Recover("dt",&dt);
    2033         if(!found)ISSMERROR(" could not find dt in inputs!");
     1952        /*recover dt: */
     1953        dt=numpar->dt;
    20341954
    20351955        /* Get node coordinates and dof list: */
    … …  
    20691989
    20701990                //Get vx, vy and their derivatives at gauss point
    2071                 GetParameterValue(&vx, &vx_list[0],gauss_l1l2l3);
    2072                 GetParameterValue(&vy, &vy_list[0],gauss_l1l2l3);
     1991                inputs->GetParameterValue(&vx,&gauss_l1l2l3[0],VxAverageEnum);
     1992                inputs->GetParameterValue(&vy,&gauss_l1l2l3[0],VyAverageEnum);
    20731993
    20741994                DL_scalar=-dt*gauss_weight*Jdettria;
    … …  
    21032023/*FUNCTION Tria::CreateKMatrixSlopeCompute {{{1*/
    21042024
    2105 void  Tria::CreateKMatrixSlopeCompute(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     2025void  Tria::CreateKMatrixSlopeCompute(Mat Kgg,int analysis_type,int sub_analysis_type){
    21062026
    21072027        /* local declarations */
    … …  
    21952115/*}}}*/
    21962116/*FUNCTION Tria::CreateKMatrixThermal {{{1*/
    2197 void  Tria::CreateKMatrixThermal(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type){
     2117void  Tria::CreateKMatrixThermal(Mat Kgg,int analysis_type,int sub_analysis_type){
    21982118
    21992119        int i,j;
    … …  
    22302150        double     tl1l2l3D[3];
    22312151        double  D_scalar;
    2232         ParameterInputs* inputs=NULL;
    22332152
    22342153        /*dynamic objects pointed to by hooks: */
    … …  
    22382157        Numpar* numpar=NULL;
    22392158
    2240         /*recover pointers: */
    2241         inputs=(ParameterInputs*)vinputs;
    2242 
    22432159        /*recover objects from hooks: */
    22442160        nodes=(Node**)hnodes.deliverp();
    … …  
    22482164
    22492165        /*recover extra inputs from users, dt: */
    2250         found=inputs->Recover("dt",&dt);
    2251         if(!found)ISSMERROR(" could not find dt in inputs!");
     2166        dt=numpar->dt;
    22522167
    22532168        /* Get node coordinates and dof list: */
    … …  
    23072222/*}}}*/
    23082223/*FUNCTION Tria::CreatePVector {{{1*/
    2309 void  Tria::CreatePVector(Vec pg,void* inputs,int analysis_type,int sub_analysis_type){
     2224void  Tria::CreatePVector(Vec pg,int analysis_type,int sub_analysis_type){
    23102225       
    23112226        /*Just branch to the correct load generator, according to the type of analysis we are carrying out: */
    23122227        if (analysis_type==ControlAnalysisEnum){
    23132228               
    2314                 CreatePVectorDiagnosticHoriz( pg,inputs,analysis_type,sub_analysis_type);
     2229                CreatePVectorDiagnosticHoriz( pg,analysis_type,sub_analysis_type);
    23152230       
    23162231        }
    … …  
    23182233                if (sub_analysis_type==HorizAnalysisEnum){
    23192234               
    2320                         CreatePVectorDiagnosticHoriz( pg,inputs,analysis_type,sub_analysis_type);
     2235                        CreatePVectorDiagnosticHoriz( pg,analysis_type,sub_analysis_type);
    23212236                }
    23222237                else ISSMERROR("%s%i%s\n","sub_analysis: ",sub_analysis_type," not supported yet");
    … …  
    23242239        else if (analysis_type==SlopecomputeAnalysisEnum){
    23252240               
    2326                 CreatePVectorSlopeCompute( pg,inputs,analysis_type,sub_analysis_type);
     2241                CreatePVectorSlopeCompute( pg,analysis_type,sub_analysis_type);
    23272242        }
    23282243        else if (analysis_type==PrognosticAnalysisEnum){
    23292244
    2330                 CreatePVectorPrognostic( pg,inputs,analysis_type,sub_analysis_type);
     2245                CreatePVectorPrognostic( pg,analysis_type,sub_analysis_type);
    23312246        }
    23322247        else if (analysis_type==Prognostic2AnalysisEnum){
    23332248
    2334                 CreatePVectorPrognostic2( pg,inputs,analysis_type,sub_analysis_type);
     2249                CreatePVectorPrognostic2( pg,analysis_type,sub_analysis_type);
    23352250        }
    23362251        else if (analysis_type==BalancedthicknessAnalysisEnum){
    23372252
    2338                 CreatePVectorBalancedthickness( pg,inputs,analysis_type,sub_analysis_type);
     2253                CreatePVectorBalancedthickness( pg,analysis_type,sub_analysis_type);
    23392254        }
    23402255        else if (analysis_type==Balancedthickness2AnalysisEnum){
    23412256
    2342                 CreatePVectorBalancedthickness2( pg,inputs,analysis_type,sub_analysis_type);
     2257                CreatePVectorBalancedthickness2( pg,analysis_type,sub_analysis_type);
    23432258        }
    23442259        else if (analysis_type==BalancedvelocitiesAnalysisEnum){
    23452260
    2346                 CreatePVectorBalancedvelocities( pg,inputs,analysis_type,sub_analysis_type);
     2261                CreatePVectorBalancedvelocities( pg,analysis_type,sub_analysis_type);
    23472262        }
    23482263        else{
    … …  
    23532268/*}}}*/
    23542269/*FUNCTION Tria::CreatePVectorBalancedthickness {{{1*/
    2355 void  Tria::CreatePVectorBalancedthickness(Vec pg ,void* vinputs,int analysis_type,int sub_analysis_type){
     2270void  Tria::CreatePVectorBalancedthickness(Vec pg ,int analysis_type,int sub_analysis_type){
    23562271
    23572272
    … …  
    23822297
    23832298        /*input parameters for structural analysis (diagnostic): */
    2384         double  accumulation_list[numgrids]={0.0};
    23852299        double  accumulation_g;
    2386         double  melting_list[numgrids]={0.0};
    23872300        double  melting_g;
    2388         double  thickness_list[numgrids]={0.0};
    2389         double  thickness_g;
    2390         int     dofs[1]={0};
    2391         int     found=0;
    23922301
    23932302        /*dynamic objects pointed to by hooks: */
    … …  
    23972306        Numpar* numpar=NULL;
    23982307
    2399         ParameterInputs* inputs=NULL;
    2400 
    2401         /*recover pointers: */
    2402         inputs=(ParameterInputs*)vinputs;
    2403 
    24042308        /*recover objects from hooks: */
    24052309        nodes=(Node**)hnodes.deliverp();
    … …  
    24072311        matice=(Matice*)hmatice.delivers();
    24082312        numpar=(Numpar*)hnumpar.delivers();
    2409 
    2410         /*recover extra inputs from users, at current convergence iteration: */
    2411         found=inputs->Recover("accumulation",&accumulation_list[0],1,dofs,numgrids,(void**)nodes);
    2412         if(!found)ISSMERROR(" could not find accumulation in inputs!");
    2413         found=inputs->Recover("melting",&melting_list[0],1,dofs,numgrids,(void**)nodes);
    2414         if(!found)ISSMERROR(" could not find melting in inputs!");
    24152313
    24162314        /* Get node coordinates and dof list: */
    … …  
    24352333                GetL(&L[0], &xyz_list[0][0], gauss_l1l2l3,numberofdofspernode);
    24362334
    2437                 /* Get accumulation, melting and thickness at gauss point */
    2438                 GetParameterValue(&accumulation_g, &accumulation_list[0],gauss_l1l2l3);
    2439                 GetParameterValue(&melting_g, &melting_list[0],gauss_l1l2l3);
     2335                /* Get accumulation, melting at gauss point */
     2336                inputs->GetParameterValue(&accumulation_g, &gauss_l1l2l3[0],AccumulationRateEnum);
     2337                inputs->GetParameterValue(&melting_g, &gauss_l1l2l3[0],MeltingRateEnum);
    24402338
    24412339                /* Add value into pe_g: */
    … …  
    24562354/*}}}*/
    24572355/*FUNCTION Tria::CreatePVectorBalancedthickness2 {{{1*/
    2458 void  Tria::CreatePVectorBalancedthickness2(Vec pg ,void* vinputs,int analysis_type,int sub_analysis_type){
     2356void  Tria::CreatePVectorBalancedthickness2(Vec pg ,int analysis_type,int sub_analysis_type){
    24592357
    24602358
    … …  
    24852383
    24862384        /*input parameters for structural analysis (diagnostic): */
    2487         double  accumulation_list[numgrids]={0.0};
    24882385        double  accumulation_g;
    2489         double  melting_list[numgrids]={0.0};
    24902386        double  melting_g;
    2491         double  dhdt_list[numgrids]={0.0};
    24922387        double  dhdt_g;
    2493         int     dofs[1]={0};
    2494         int     found=0;
    24952388
    24962389        /*dynamic objects pointed to by hooks: */
    … …  
    25002393        Numpar* numpar=NULL;
    25012394
    2502         ParameterInputs* inputs=NULL;
    2503 
    2504         /*recover pointers: */
    2505         inputs=(ParameterInputs*)vinputs;
    2506 
    25072395        /*recover objects from hooks: */
    25082396        nodes=(Node**)   hnodes.deliverp();
    … …  
    25102398        matice=(Matice*)hmatice.delivers();
    25112399        numpar=(Numpar*)hnumpar.delivers();
    2512 
    2513         /*recover extra inputs from users, at current convergence iteration: */
    2514         found=inputs->Recover("accumulation",&accumulation_list[0],1,dofs,numgrids,(void**)nodes);
    2515         if(!found)ISSMERROR(" could not find accumulation in inputs!");
    2516         found=inputs->Recover("melting",&melting_list[0],1,dofs,numgrids,(void**)nodes);
    2517         if(!found)ISSMERROR(" could not find melting in inputs!");
    2518         found=inputs->Recover("dhdt",&dhdt_list[0],1,dofs,numgrids,(void**)nodes);
    2519         if(!found)ISSMERROR(" could not find dhdt in inputs!");
    25202400
    25212401        /* Get node coordinates and dof list: */
    … …  
    25412421
    25422422                /* Get accumulation, melting and thickness at gauss point */
    2543                 GetParameterValue(&accumulation_g, &accumulation_list[0],gauss_l1l2l3);
    2544                 GetParameterValue(&melting_g, &melting_list[0],gauss_l1l2l3);
    2545                 GetParameterValue(&dhdt_g, &dhdt_list[0],gauss_l1l2l3);
     2423                inputs->GetParameterValue(&accumulation_g, &gauss_l1l2l3[0],AccumulationRateEnum);
     2424                inputs->GetParameterValue(&melting_g, &gauss_l1l2l3[0],MeltingRateEnum);
     2425                inputs->GetParameterValue(&dhdt_g, &gauss_l1l2l3[0],DhDtEnum);
    25462426
    25472427                /* Add value into pe_g: */
    … …  
    25622442/*}}}*/
    25632443/*FUNCTION Tria::CreatePVectorBalancedvelocities {{{1*/
    2564 void  Tria::CreatePVectorBalancedvelocities(Vec pg ,void* vinputs,int analysis_type,int sub_analysis_type){
     2444void  Tria::CreatePVectorBalancedvelocities(Vec pg ,int analysis_type,int sub_analysis_type){
    25652445
    25662446
    … …  
    25912471
    25922472        /*input parameters for structural analysis (diagnostic): */
    2593         double  accumulation_list[numgrids]={0.0};
    25942473        double  accumulation_g;
    2595         double  melting_list[numgrids]={0.0};
    25962474        double  melting_g;
    2597         int     dofs[1]={0};
    2598         int     found=0;
    25992475
    26002476        /*dynamic objects pointed to by hooks: */
    … …  
    26042480        Numpar* numpar=NULL;
    26052481
    2606         ParameterInputs* inputs=NULL;
    2607 
    2608         /*recover pointers: */
    2609         inputs=(ParameterInputs*)vinputs;
    2610 
    26112482        /*recover objects from hooks: */
    26122483        nodes=(Node**)hnodes.deliverp();
    … …  
    26202491        matice=(Matice*)hmatice.delivers();
    26212492        numpar=(Numpar*)hnumpar.delivers();
    2622 
    2623         /*recover extra inputs from users, at current convergence iteration: */
    2624         found=inputs->Recover("accumulation",&accumulation_list[0],1,dofs,numgrids,(void**)nodes);
    2625         if(!found)ISSMERROR(" could not find accumulation in inputs!");
    2626         found=inputs->Recover("melting",&melting_list[0],1,dofs,numgrids,(void**)nodes);
    2627         if(!found)ISSMERROR(" could not find melting in inputs!");
    26282493
    26292494        /* Get node coordinates and dof list: */
    … …  
    26492514
    26502515                /* Get accumulation, melting at gauss point */
    2651                 GetParameterValue(&accumulation_g, &accumulation_list[0],gauss_l1l2l3);
    2652                 GetParameterValue(&melting_g, &melting_list[0],gauss_l1l2l3);
     2516                inputs->GetParameterValue(&accumulation_g, &gauss_l1l2l3[0],AccumulationRateEnum);
     2517                inputs->GetParameterValue(&melting_g, &gauss_l1l2l3[0],MeltingRateEnum);
    26532518
    26542519                /* Add value into pe_g: */
    … …  
    26692534/*}}}*/
    26702535/*FUNCTION Tria::CreatePVectorDiagnosticBaseVert {{{1*/
    2671 void  Tria::CreatePVectorDiagnosticBaseVert(Vec pg,void* vinputs,int analysis_type,int sub_analysis_type){
     2536void  Tria::CreatePVectorDiagnosticBaseVert(Vec pg,int analysis_type,int sub_analysis_type){
    26722537
    26732538        int             i,j;
    … …  
    27042569
    27052570        /*input parameters for structural analysis (diagnostic): */
    2706         double* velocity_param=NULL;
    2707         double  vx_list[numgrids]={0,0,0};
    2708         double  vy_list[numgrids]={0,0,0};
    27092571        double  vx,vy;
    27102572        double  meltingvalue;
    27112573        double  slope[2];
    27122574        double  dbdx,dbdy;
    2713         int     dofs1[1]={0};
    2714         int     dofs2[1]={1};
    27152575
    27162576        /*dynamic objects pointed to by hooks: */
    … …  
    27202580        Numpar* numpar=NULL;
    27212581
    2722         ParameterInputs* inputs=NULL;
    2723 
    2724         /*recover pointers: */
    2725         inputs=(ParameterInputs*)vinputs;
    2726 
    27272582        /*recover objects from hooks: */
    27282583        nodes=(Node**)hnodes.deliverp();
    … …  
    27312586        numpar=(Numpar*)hnumpar.delivers();
    27322587
    2733         /* recover input parameters: */
    2734         if(!inputs->Recover("velocity",&vx_list[0],1,dofs1,numgrids,(void**)nodes))ISSMERROR(" cannot compute vertical velocity without horizontal velocity");
    2735         inputs->Recover("velocity",&vy_list[0],1,dofs2,numgrids,(void**)nodes);
    2736 
    27372588        /* Get node coordinates and dof list: */
    27382589        GetVerticesCoordinates(&xyz_list[0][0], nodes, numgrids);
    … …  
    27562607
    27572608                /*Get melting at gaussian point: */
    2758                 GetParameterValue(&meltingvalue, &this->properties.melting[0],gauss_l1l2l3);
     2609                inputs->GetParameterValue(&meltingvalue, &gauss_l1l2l3[0],MeltingRateEnum);
    27592610
    27602611                /*Get velocity at gaussian point: */
    2761                 GetParameterValue(&vx, &vx_list[0],gauss_l1l2l3);
    2762                 GetParameterValue(&vy, &vy_list[0],gauss_l1l2l3);
     2612                inputs->GetParameterValue(&vx, &gauss_l1l2l3[0],VxEnum);
     2613                inputs->GetParameterValue(&vy, &gauss_l1l2l3[0],VyEnum);
    27632614
    27642615                /*Get bed slope: */
    2765                 GetParameterDerivativeValue(&slope[0], &this->properties.b[0],&xyz_list[0][0], gauss_l1l2l3);
     2616                inputs->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],&gauss_l1l2l3[0],BedEnum);
    27662617                dbdx=slope[0];
    27672618                dbdy=slope[1];
    … …  
    27962647/*}}}*/
    27972648/*FUNCTION Tria::CreatePVectorDiagnosticHoriz {{{1*/
    2798 void Tria::CreatePVectorDiagnosticHoriz( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type){
     2649void Tria::CreatePVectorDiagnosticHoriz( Vec pg, int analysis_type,int sub_analysis_type){
    27992650
    28002651        int             i,j;
    … …  
    28412692        Numpar* numpar=NULL;
    28422693
    2843         ParameterInputs* inputs=NULL;
     2694        /*inputs: */
     2695        bool onwater;
     2696        int  drag_type;
     2697
     2698        /*retrieve inputs :*/
     2699        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     2700        inputs->GetParameterValue(&drag_type,DragTypeEnum);
    28442701
    28452702        /*First, if we are on water, return empty vector: */
    2846         if(this->properties.onwater)return;
    2847 
    2848         /*recover pointers: */
    2849         inputs=(ParameterInputs*)vinputs;
     2703        if(onwater)return;
    28502704
    28512705        /*recover objects from hooks: */
    … …  
    28752729
    28762730                /*Compute thickness at gaussian point: */
    2877                 GetParameterValue(&thickness, &this->properties.h[0],gauss_l1l2l3);
    2878        
    2879                 GetParameterDerivativeValue(&slope[0], &this->properties.s[0],&xyz_list[0][0], gauss_l1l2l3);
     2731                inputs->GetParameterValue(&thickness, &gauss_l1l2l3[0],ThicknessEnum);
     2732                inputs->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],&gauss_l1l2l3[0],SurfaceEnum);
    28802733               
    28812734                /*In case we have plastic basal drag, compute plastic stress at gaussian point from k1, k2 and k3 fields in the
    28822735                 * element itself: */
    2883                 if(this->properties.friction_type==1){
    2884                         GetParameterValue(&plastic_stress, &this->properties.k[0],gauss_l1l2l3);
     2736                if(drag_type==1){
     2737                        inputs->GetParameterValue(&plastic_stress, &gauss_l1l2l3[0],DragCoefficientEnum);
    28852738                }
    28862739
    … …  
    28952748
    28962749                /*Build pe_g_gaussian vector: */
    2897                 if(this->properties.friction_type==1){
     2750                if(drag_type==1){
    28982751                        for (i=0;i<numgrids;i++){
    28992752                                for (j=0;j<NDOF2;j++){
    … …  
    29272780/*}}}*/
    29282781/*FUNCTION Tria::CreatePVectorPrognostic {{{1*/
    2929 void  Tria::CreatePVectorPrognostic(Vec pg ,void* vinputs,int analysis_type,int sub_analysis_type){
     2782void  Tria::CreatePVectorPrognostic(Vec pg ,int analysis_type,int sub_analysis_type){
    29302783
    29312784
    … …  
    29562809
    29572810        /*input parameters for structural analysis (diagnostic): */
    2958         double  accumulation_list[numgrids]={0.0};
    29592811        double  accumulation_g;
    2960         double  melting_list[numgrids]={0.0};
    29612812        double  melting_g;
    2962         double  thickness_list[numgrids]={0.0};
    29632813        double  thickness_g;
    29642814        double  dt;
    2965         int     dofs[1]={0};
    2966         int     found=0;
    29672815
    29682816        /*dynamic objects pointed to by hooks: */
    … …  
    29722820        Numpar* numpar=NULL;
    29732821
    2974         ParameterInputs* inputs=NULL;
    2975 
    2976         /*recover pointers: */
    2977         inputs=(ParameterInputs*)vinputs;
    2978 
    29792822        /*recover objects from hooks: */
    29802823        nodes=(Node**)hnodes.deliverp();
    … …  
    29832826        numpar=(Numpar*)hnumpar.delivers();
    29842827
    2985         /*recover extra inputs from users, at current convergence iteration: */
    2986         found=inputs->Recover("accumulation",&accumulation_list[0],1,dofs,numgrids,(void**)nodes);
    2987         if(!found)ISSMERROR(" could not find accumulation in inputs!");
    2988         found=inputs->Recover("melting",&melting_list[0],1,dofs,numgrids,(void**)nodes);
    2989         if(!found)ISSMERROR(" could not find melting in inputs!");
    2990         found=inputs->Recover("thickness",&thickness_list[0],1,dofs,numgrids,(void**)nodes);
    2991         if(!found)ISSMERROR(" could not find thickness in inputs!");
    2992         found=inputs->Recover("dt",&dt);
    2993         if(!found)ISSMERROR(" could not find dt in inputs!");
     2828        /*recover dt: */
     2829        dt=numpar->dt;
    29942830
    29952831        /* Get node coordinates and dof list: */
    … …  
    30152851
    30162852                /* Get accumulation, melting and thickness at gauss point */
    3017                 GetParameterValue(&accumulation_g, &accumulation_list[0],gauss_l1l2l3);
    3018                 GetParameterValue(&melting_g, &melting_list[0],gauss_l1l2l3);
    3019                 GetParameterValue(&thickness_g, &thickness_list[0],gauss_l1l2l3);
     2853                inputs->GetParameterValue(&accumulation_g, &gauss_l1l2l3[0],AccumulationRateEnum);
     2854                inputs->GetParameterValue(&melting_g, &gauss_l1l2l3[0],MeltingRateEnum);
     2855                inputs->GetParameterValue(&thickness_g, &gauss_l1l2l3[0],ThicknessEnum);
    30202856
    30212857                /* Add value into pe_g: */
    … …  
    30362872/*}}}*/
    30372873/*FUNCTION Tria::CreatePVectorPrognostic2 {{{1*/
    3038 void  Tria::CreatePVectorPrognostic2(Vec pg ,void* vinputs,int analysis_type,int sub_analysis_type){
     2874void  Tria::CreatePVectorPrognostic2(Vec pg ,int analysis_type,int sub_analysis_type){
    30392875
    30402876        /* local declarations */
    … …  
    30642900
    30652901        /*input parameters for structural analysis (diagnostic): */
    3066         double  accumulation_list[numgrids]={0.0};
    30672902        double  accumulation_g;
    3068         double  melting_list[numgrids]={0.0};
    30692903        double  melting_g;
    3070         double  thickness_list[numgrids]={0.0};
    30712904        double  thickness_g;
    30722905        double  dt;
    3073         int     dofs[1]={0};
    3074         int     found=0;
    30752906
    30762907        /*dynamic objects pointed to by hooks: */
    … …  
    30802911        Numpar* numpar=NULL;
    30812912
    3082         ParameterInputs* inputs=NULL;
    3083 
    3084         /*recover pointers: */
    3085         inputs=(ParameterInputs*)vinputs;
    3086 
    30872913        /*recover objects from hooks: */
    30882914        nodes=(Node**)hnodes.deliverp();
    … …  
    30912917        numpar=(Numpar*)hnumpar.delivers();
    30922918
    3093         /*recover extra inputs from users, at current convergence iteration: */
    3094         found=inputs->Recover("accumulation",&accumulation_list[0],1,dofs,numgrids,(void**)nodes);
    3095         if(!found)ISSMERROR(" could not find accumulation in inputs!");
    3096         found=inputs->Recover("melting",&melting_list[0],1,dofs,numgrids,(void**)nodes);
    3097         if(!found)ISSMERROR(" could not find melting in inputs!");
    3098         found=inputs->Recover("thickness",&thickness_list[0],1,dofs,numgrids,(void**)nodes);
    3099         if(!found)ISSMERROR(" could not find thickness in inputs!");
    3100         found=inputs->Recover("dt",&dt);
    3101         if(!found)ISSMERROR(" could not find dt in inputs!");
     2919        /*recover dt: */
     2920        dt=numpar->dt;
     2921
    31022922
    31032923        /* Get node coordinates and dof list: */
    … …  
    31232943
    31242944                /* Get accumulation, melting and thickness at gauss point */
    3125                 GetParameterValue(&accumulation_g, &accumulation_list[0],gauss_l1l2l3);
    3126                 GetParameterValue(&melting_g, &melting_list[0],gauss_l1l2l3);
    3127                 GetParameterValue(&thickness_g, &thickness_list[0],gauss_l1l2l3);
     2945                inputs->GetParameterValue(&accumulation_g, &gauss_l1l2l3[0],AccumulationRateEnum);
     2946                inputs->GetParameterValue(&melting_g, &gauss_l1l2l3[0],MeltingRateEnum);
     2947                inputs->GetParameterValue(&thickness_g, &gauss_l1l2l3[0],ThicknessEnum);
    31282948
    31292949                /* Add value into pe_g: */
    … …  
    31452965/*FUNCTION Tria::CreatePVectorSlopeCompute {{{1*/
    31462966
    3147 void Tria::CreatePVectorSlopeCompute( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type){
     2967void Tria::CreatePVectorSlopeCompute( Vec pg, int analysis_type,int sub_analysis_type){
    31482968
    31492969        int             i,j;
    … …  
    31752995        double  pe_g[numdof];
    31762996        double  pe_g_gaussian[numdof];
    3177         double  param[3];
    31782997        double  slope[2];
    31792998
    … …  
    31973016        for(i=0;i<numdof;i++) pe_g[i]=0.0;
    31983017
    3199         if ( (sub_analysis_type==SurfaceXAnalysisEnum) || (sub_analysis_type==SurfaceYAnalysisEnum)){
    3200                 for(i=0;i<numdof;i++) param[i]=this->properties.s[i];
    3201         }
    3202         if ( (sub_analysis_type==BedXAnalysisEnum) || (sub_analysis_type==BedYAnalysisEnum)){
    3203                 for(i=0;i<numdof;i++) param[i]=this->properties.b[i];
    3204         }
    32053018
    32063019        /* Get gaussian points and weights: */
    … …  
    32163029                gauss_l1l2l3[2]=*(third_gauss_area_coord+ig);
    32173030
    3218                 GetParameterDerivativeValue(&slope[0], &param[0],&xyz_list[0][0], gauss_l1l2l3);
    3219                
     3031                if ( (sub_analysis_type==SurfaceXAnalysisEnum) || (sub_analysis_type==SurfaceYAnalysisEnum)){
     3032                        inputs->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],&gauss_l1l2l3[0],SurfaceEnum);
     3033                }
     3034                if ( (sub_analysis_type==BedXAnalysisEnum) || (sub_analysis_type==BedYAnalysisEnum)){
     3035                        inputs->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],&gauss_l1l2l3[0],BedEnum);
     3036                }
     3037
    32203038                /* Get Jacobian determinant: */
    32213039                GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0],gauss_l1l2l3);
    … …  
    32493067/*}}}*/
    32503068/*FUNCTION Tria::CreatePVectorThermalShelf {{{1*/
    3251 void Tria::CreatePVectorThermalShelf( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type){
     3069void Tria::CreatePVectorThermalShelf( Vec pg, int analysis_type,int sub_analysis_type){
    32523070
    32533071        int i,found;
    … …  
    32703088        /*inputs: */
    32713089        double dt;
    3272         double pressure_list[3];
    32733090        double pressure;
    32743091
    … …  
    32973114        Numpar* numpar=NULL;
    32983115
    3299         ParameterInputs* inputs=NULL;
    3300 
    3301         /*recover pointers: */
    3302         inputs=(ParameterInputs*)vinputs;
    3303        
    33043116        /*recover objects from hooks: */
    33053117        nodes=(Node**)hnodes.deliverp();
    … …  
    33203132        beta=matpar->GetBeta();
    33213133        meltingpoint=matpar->GetMeltingPoint();
    3322 
    3323 
    3324         /*recover extra inputs from users, dt and velocity: */
    3325         found=inputs->Recover("dt",&dt);
    3326         if(!found)ISSMERROR(" could not find dt in inputs!");
    3327         found=inputs->Recover("pressure",&pressure_list[0],1,dofs1,numgrids,(void**)nodes);
    3328         if(!found)ISSMERROR(" could not find pressure in inputs!");
     3134        dt=numpar->dt;
     3135
    33293136
    33303137        /* Ice/ocean heat exchange flux on ice shelf base */
    … …  
    33463153
    33473154                /*Get geothermal flux and basal friction */
    3348                 GetParameterValue(&pressure,&pressure_list[0],gauss_coord);
     3155                inputs->GetParameterValue(&pressure, &gauss_coord[0],PressureEnum);
    33493156                t_pmp=meltingpoint-beta*pressure;
    33503157
    … …  
    33723179/*}}}*/
    33733180/*FUNCTION Tria::CreatePVectorThermalSheet {{{1*/
    3374 void Tria::CreatePVectorThermalSheet( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type){
     3181void Tria::CreatePVectorThermalSheet( Vec pg, int analysis_type,int sub_analysis_type){
    33753182
    33763183        int i,found;
    … …  
    33823189        int        numberofdofspernode;
    33833190        double       xyz_list[numgrids][3];
    3384         double     vxvyvz_list[numgrids][3];
    3385         double     vx_list[numgrids];
    3386         double     vy_list[numgrids];
    33873191
    33883192        double rho_ice;
    … …  
    33983202        double geothermalflux_value;
    33993203
     3204        double  vx_list[numgrids];
     3205        double  vy_list[numgrids];
     3206        double  thickness_list[numgrids];
     3207        double  bed_list[numgrids];
     3208        double  dragcoefficient_list[numgrids];
     3209        double  drag_p,drag_q;
     3210        int     drag_type;
     3211
    34003212        /* gaussian points: */
    34013213        int     num_area_gauss,ig;
    … …  
    34063218        double  gauss_weight;
    34073219        double  gauss_coord[3];
    3408         int     dofs1[1]={0};
     3220        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
    34093221
    34103222        /*matrices: */
    … …  
    34143226        double  scalar;
    34153227
    3416         int     dofs[3]={0,1,2};
    3417 
    3418         ParameterInputs* inputs=NULL;
    34193228
    34203229        /*dynamic objects pointed to by hooks: */
    … …  
    34243233        Numpar* numpar=NULL;
    34253234
    3426         /*recover pointers: */
    3427         inputs=(ParameterInputs*)vinputs;
    3428 
    34293235        /*recover objects from hooks: */
    34303236        nodes=(Node**)hnodes.deliverp();
    … …  
    34323238        matice=(Matice*)hmatice.delivers();
    34333239        numpar=(Numpar*)hnumpar.delivers();
     3240
     3241        /*retrieve inputs :*/
     3242        inputs->GetParameterValue(&drag_type,DragTypeEnum);
    34343243       
    34353244        /* Get node coordinates and dof list: */
    … …  
    34403249        rho_ice=matpar->GetRhoIce();
    34413250        heatcapacity=matpar->GetHeatCapacity();
    3442 
    3443 
    3444         /*recover extra inputs from users, dt and velocity: */
    3445         found=inputs->Recover("dt",&dt);
    3446         if(!found)ISSMERROR(" could not find dt in inputs!");
    3447        
    3448         found=inputs->Recover("velocity",&vxvyvz_list[0][0],3,dofs,numgrids,(void**)nodes);
    3449         if(!found)ISSMERROR(" could not find velocity in inputs!");
    3450 
    3451         for(i=0;i<numgrids;i++){
    3452                 vx_list[i]=vxvyvz_list[i][0];
    3453                 vy_list[i]=vxvyvz_list[i][1];
    3454         }       
     3251        dt=numpar->dt;
     3252
     3253
     3254        /*Recover inputs: */
     3255        inputs->GetParameterValues(&vx_list[0],&gaussgrids[0][0],3,VxAverageEnum);
     3256        inputs->GetParameterValues(&vy_list[0],&gaussgrids[0][0],3,VyAverageEnum);
     3257        inputs->GetParameterValues(&dragcoefficient_list[0],&gaussgrids[0][0],3,DragCoefficientEnum);
     3258        inputs->GetParameterValues(&bed_list[0],&gaussgrids[0][0],3,BedEnum);
     3259        inputs->GetParameterValues(&thickness_list[0],&gaussgrids[0][0],3,ThicknessEnum);
     3260        inputs->GetParameterValue(&drag_p,DragPEnum);
     3261        inputs->GetParameterValue(&drag_q,DragQEnum);
    34553262
    34563263        /*Build alpha2_list used by drag stiffness matrix*/
    … …  
    34583265       
    34593266        /*Initialize all fields: */
    3460         if (this->properties.friction_type!=2)ISSMERROR(" non-viscous friction not supported yet!");
     3267        if (drag_type!=2)ISSMERROR(" non-viscous friction not supported yet!");
    34613268       
    34623269        friction->element_type=(char*)xmalloc((strlen("3d")+1)*sizeof(char));
    34633270        strcpy(friction->element_type,"3d");
    3464        
     3271
    34653272        friction->gravity=matpar->GetG();
    34663273        friction->rho_ice=matpar->GetRhoIce();
    34673274        friction->rho_water=matpar->GetRhoWater();
    3468         friction->K=&this->properties.k[0];
    3469         friction->bed=&this->properties.b[0];
    3470         friction->thickness=&this->properties.h[0];
    3471         friction->velocities=&vxvyvz_list[0][0];
    3472         friction->p=this->properties.p;
    3473         friction->q=this->properties.q;
     3275        friction->K=&dragcoefficient_list[0];
     3276        friction->bed=&bed_list[0];
     3277        friction->thickness=&thickness_list[0];
     3278        friction->vx=&vx_list[0];
     3279        friction->vy=&vy_list[0];
     3280        friction->p=drag_p;
     3281        friction->q=drag_q;
    34743282
    34753283        /*Compute alpha2_list: */
    … …  
    35013309
    35023310                /*Get geothermal flux and basal friction */
    3503                 GetParameterValue(&geothermalflux_value,&this->properties.geothermalflux[0],gauss_coord);
     3311                inputs->GetParameterValue(&geothermalflux_value, &gauss_coord[0],GeothermalFluxEnum);
    35043312                GetParameterValue(&basalfriction,&basalfriction_list[0],gauss_coord);
    35053313
    … …  
    35273335/*}}}*/
    35283336/*FUNCTION Tria::Du {{{1*/
    3529 void Tria::Du(Vec du_g,void* vinputs,int analysis_type,int sub_analysis_type){
     3337void Tria::Du(Vec du_g,int analysis_type,int sub_analysis_type){
    35303338
    35313339        int i;
    … …  
    35383346        int          doflist[numdof];
    35393347        int          numberofdofspernode;
    3540         int          dofs2[2]={0,1};
    3541         int          dofs1[1]={0};
     3348        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
    35423349
    35433350        /* grid data: */
    3544         double vxvy_list[numgrids][2];
    35453351        double vx_list[numgrids];
    35463352        double vy_list[numgrids];
    3547         double obs_vxvy_list[numgrids][2];
    35483353        double obs_vx_list[numgrids];
    35493354        double obs_vy_list[numgrids];
    … …  
    35663371
    35673372        /*element vector : */
    3568         double  due_g[numdof];
     3373        double  due_g[numdof]={0,0,0,0,0,0};
    35693374        double  due_g_gaussian[numdof];
    35703375
    … …  
    35883393        Numpar* numpar=NULL;
    35893394
    3590         ParameterInputs* inputs=NULL;
    3591 
    3592         /*recover pointers: */
    3593         inputs=(ParameterInputs*)vinputs;
    3594 
    35953395        /*recover objects from hooks: */
    35963396        nodes=(Node**)hnodes.deliverp();
    … …  
    36033403        GetDofList(&doflist[0],&numberofdofspernode);
    36043404
    3605         /* Set due_g to 0: */
    3606         for(i=0;i<numdof;i++) due_g[i]=0.0;
    3607 
    36083405        /* Recover input data: */
    3609         if(!inputs->Recover("fit",&fit)) ISSMERROR(" missing fit input parameter");
    3610         if(!inputs->Recover("velocity_obs",&obs_vxvy_list[0][0],2,dofs2,numgrids,(void**)nodes)){
    3611                 ISSMERROR("missing velocity_obs input parameter");
    3612         }
    3613         if(fit==3 && !inputs->Recover("surfacearea",&S)){
    3614                 ISSMERROR("missing surface area input parameter");
    3615         }
    3616         if(!inputs->Recover("velocity",&vxvy_list[0][0],2,dofs2,numgrids,(void**)nodes)){
    3617                 ISSMERROR("missing velocity input parameter");
    3618         }
    3619         if(!inputs->Recover("weights",&weights_list[0],1,dofs1,numgrids,(void**)nodes)){
    3620                 ISSMERROR("missing weights input parameter");
    3621         }
    3622 
    3623         for(i=0;i<numgrids;i++){
    3624                 obs_vx_list[i]=obs_vxvy_list[i][0];
    3625                 obs_vy_list[i]=obs_vxvy_list[i][1];
    3626                 vx_list[i]=vxvy_list[i][0];
    3627                 vy_list[i]=vxvy_list[i][1];
     3406        inputs->GetParameterValues(&obs_vx_list[0],&gaussgrids[0][0],3,VxObsEnum);
     3407        inputs->GetParameterValues(&obs_vy_list[0],&gaussgrids[0][0],3,VyObsEnum);
     3408
     3409        inputs->GetParameterValues(&vx_list[0],&gaussgrids[0][0],3,VxEnum);
     3410        inputs->GetParameterValues(&vy_list[0],&gaussgrids[0][0],3,VyEnum);
     3411       
     3412        inputs->GetParameterValues(&weights_list[0],&gaussgrids[0][0],3,WeightsEnum);
     3413       
     3414        inputs->GetParameterValue(&fit,FitEnum);
     3415        if(fit==3){
     3416                inputs->GetParameterValue(&S,SurfaceAreaEnum);
    36283417        }
    36293418
    … …  
    39273716/*}}}*/
    39283717/*FUNCTION Tria::GetBedList {{{1*/
    3929 void  Tria::GetBedList(double* bed_list){
    3930 
    3931         int i;
    3932         for(i=0;i<3;i++)bed_list[i]=this->properties.b[i];
     3718void  Tria::GetBedList(double* bedlist){
     3719       
     3720        const int numgrids=3;
     3721        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
     3722       
     3723        inputs->GetParameterValues(bedlist,&gaussgrids[0][0],3,BedEnum);
    39333724
    39343725}
    … …  
    42964087/*}}}*/
    42974088/*FUNCTION Tria::GetOnBed {{{1*/
    4298 int Tria::GetOnBed(){
    4299         return this->properties.onbed;
     4089bool Tria::GetOnBed(){
     4090       
     4091        bool onbed;
     4092        inputs->GetParameterValue(&onbed,ElementOnBedEnum);
     4093
     4094        return onbed;
    43004095}
    43014096/*}}}*/
    … …  
    43444139/*}}}*/
    43454140/*FUNCTION Tria::GetShelf {{{1*/
    4346 int   Tria::GetShelf(){
    4347         return this->properties.shelf;
     4141bool   Tria::GetShelf(){
     4142        /*inputs: */
     4143        bool shelf;
     4144
     4145        /*retrieve inputs :*/
     4146        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
     4147
     4148        return shelf;
    43484149}
    43494150/*}}}*/
    … …  
    43704171void Tria::GetThicknessList(double* thickness_list){
    43714172
    4372         int i;
    4373         for(i=0;i<3;i++)thickness_list[i]=this->properties.h[i];
     4173        const int numgrids=3;
     4174        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
     4175        inputs->GetParameterValues(thickness_list,&gaussgrids[0][0],3,ThicknessEnum);
     4176
    43744177}
    43754178/*}}}*/
    43764179/*FUNCTION Tria::Gradj {{{1*/
    4377 void  Tria::Gradj(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type,char* control_type){
     4180void  Tria::Gradj(Vec grad_g,int analysis_type,int sub_analysis_type,char* control_type){
     4181
     4182        /*inputs: */
     4183        bool onwater;
     4184
     4185        /*retrieve inputs :*/
     4186        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
    43784187
    43794188        /*If on water, grad = 0: */
    4380         if(this->properties.onwater)return;
     4189        if(onwater)return;
    43814190
    43824191        if (strcmp(control_type,"drag")==0){
    4383                 GradjDrag( grad_g,inputs,analysis_type,sub_analysis_type);
     4192                GradjDrag( grad_g,analysis_type,sub_analysis_type);
    43844193        }
    43854194        else if (strcmp(control_type,"B")==0){
    4386                 GradjB( grad_g,inputs,analysis_type,sub_analysis_type);
     4195                GradjB( grad_g,analysis_type,sub_analysis_type);
    43874196        }
    43884197        else ISSMERROR("%s%s","control type not supported yet: ",control_type);
    … …  
    43904199/*}}}*/
    43914200/*FUNCTION Tria::GradjB{{{1*/
    4392 void  Tria::GradjB(Vec grad_g,void* vinputs,int analysis_type,int sub_analysis_type){
     4201void  Tria::GradjB(Vec grad_g,int analysis_type,int sub_analysis_type){
    43934202
    43944203        int i;
    … …  
    44044213
    44054214        /* grid data: */
    4406         double vx_list[numgrids];
    4407         double vy_list[numgrids];
    4408         double vxvy_list[numgrids][2];
    4409         double adjx_list[numgrids];
    4410         double adjy_list[numgrids];
    4411         double adjxadjy_list[numgrids][2];
    44124215        double B[numgrids];
    44134216
    … …  
    44564259        Numpar* numpar=NULL;
    44574260
    4458         ParameterInputs* inputs=NULL;
    4459 
    4460         /*recover pointers: */
    4461         inputs=(ParameterInputs*)vinputs;
    4462 
    44634261        /*recover objects from hooks: */
    44644262        nodes=(Node**)hnodes.deliverp();
    … …  
    44734271        /* Set grade_g to 0: */
    44744272        for(i=0;i<numgrids;i++) grade_g[i]=0.0;
    4475 
    4476         /* recover input parameters: */
    4477         inputs->Recover("thickness",&this->properties.h[0],1,dofs,numgrids,(void**)nodes);
    4478         if(!inputs->Recover("velocity",&vxvy_list[0][0],2,dofs2,numgrids,(void**)nodes)){
    4479                 ISSMERROR("missing velocity input parameter");
    4480         }
    4481         if(!inputs->Recover("adjoint",&adjxadjy_list[0][0],2,dofs2,numgrids,(void**)nodes)){
    4482                 ISSMERROR("missing adjoint input parameter");
    4483         }
    4484         if(!inputs->Recover("B",&B[0],1,dofs1,numgrids,(void**)nodes)){
    4485                 ISSMERROR("parameter B not found in input");
    4486         }
    4487 
    4488         /*Initialize parameter lists: */
    4489         for(i=0;i<numgrids;i++){
    4490                 vx_list[i]=vxvy_list[i][0];
    4491                 vy_list[i]=vxvy_list[i][1];
    4492                 adjx_list[i]=adjxadjy_list[i][0];
    4493                 adjy_list[i]=adjxadjy_list[i][1];
    4494         }
    44954273
    44964274        /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
    … …  
    45064284
    45074285                /*Get thickness: */
    4508                 GetParameterValue(&thickness, &this->properties.h[0],gauss_l1l2l3);
     4286                inputs->GetParameterValue(&thickness, gauss_l1l2l3,ThicknessEnum);
    45094287
    45104288                /*Get strain rate, if velocity has been supplied: */
    4511                 GetStrainRate(&epsilon[0],&vxvy_list[0][0],&xyz_list[0][0],gauss_l1l2l3);
     4289                inputs->GetStrainRate(&epsilon[0],&xyz_list[0][0],gauss_l1l2l3,VxEnum,VyEnum);
    45124290
    45134291                /*Get viscosity complement: */
    … …  
    45154293
    45164294                /*Get dvx, dvy, dadjx and dadjx: */
    4517                 GetParameterDerivativeValue(&dvx[0], &vx_list[0],&xyz_list[0][0], gauss_l1l2l3);
    4518                 GetParameterDerivativeValue(&dvy[0], &vy_list[0],&xyz_list[0][0], gauss_l1l2l3);
    4519                 GetParameterDerivativeValue(&dadjx[0], &adjx_list[0],&xyz_list[0][0], gauss_l1l2l3);
    4520                 GetParameterDerivativeValue(&dadjy[0], &adjy_list[0],&xyz_list[0][0], gauss_l1l2l3);
     4295                inputs->GetParameterDerivativeValue(&dvx[0],&xyz_list[0][0],&gauss_l1l2l3[0],VxEnum);
     4296                inputs->GetParameterDerivativeValue(&dvy[0],&xyz_list[0][0],&gauss_l1l2l3[0],VyEnum);
     4297                inputs->GetParameterDerivativeValue(&dadjx[0],&xyz_list[0][0],&gauss_l1l2l3[0],AdjointxEnum);
     4298                inputs->GetParameterDerivativeValue(&dadjy[0],&xyz_list[0][0],&gauss_l1l2l3[0],AdjointyEnum);
    45214299
    45224300                /* Get Jacobian determinant: */
    … …  
    45304308
    45314309                /*Get B derivative: dB/dx */
    4532                 GetParameterDerivativeValue(&dB[0], &B[0],&xyz_list[0][0], gauss_l1l2l3);
    4533                 GetParameterValue(&B_gauss, &B[0],gauss_l1l2l3);
     4310                inputs->GetParameterDerivativeValue(&dB[0],&xyz_list[0][0],&gauss_l1l2l3[0],RheologyBEnum);
     4311                inputs->GetParameterValue(&B_gauss, gauss_l1l2l3,RheologyBEnum);
    45344312
    45354313                /*Build gradje_g_gaussian vector (actually -dJ/dB): */
    … …  
    45684346/*}}}*/
    45694347/*FUNCTION Tria::GradjDrag {{{1*/
    4570 void  Tria::GradjDrag(Vec grad_g,void* vinputs,int analysis_type,int sub_analysis_type){
     4348void  Tria::GradjDrag(Vec grad_g,int analysis_type,int sub_analysis_type){
    45714349
    45724350
    … …  
    45844362        double vx_list[numgrids];
    45854363        double vy_list[numgrids];
    4586         double vxvy_list[numgrids][2];
    45874364        double adjx_list[numgrids];
    45884365        double adjy_list[numgrids];
    4589         double adjxadjy_list[numgrids][2];
     4366        double thickness_list[numgrids];
     4367        double bed_list[numgrids];
     4368        double dragcoefficient_list[numgrids];
     4369        double drag_p;
     4370        double drag_q;
     4371        int    drag_type;
    45904372
    45914373        double drag;
    4592         int    dofs1[1]={0};
    4593         int    dofs2[2]={0,1};
    45944374
    45954375        /* gaussian points: */
    … …  
    46014381        double  gauss_weight;
    46024382        double  gauss_l1l2l3[3];
     4383        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
    46034384
    46044385        /* parameters: */
    … …  
    46094390       
    46104391        /*drag: */
    4611         double  pcoeff,qcoeff;
    46124392        double alpha_complement_list[numgrids];
    46134393        double alpha_complement;
    46144394
    46154395        /*element vector at the gaussian points: */
    4616         double  grade_g[numgrids];
     4396        double  grade_g[numgrids]={0,0,0};
    46174397        double  grade_g_gaussian[numgrids];
    46184398
    … …  
    46254405        /* strain rate: */
    46264406        double epsilon[3]; /* epsilon=[exx,eyy,exy];*/
     4407
     4408        /*inputs: */
     4409        bool shelf;
    46274410
    46284411        /*dynamic objects pointed to by hooks: */
    … …  
    46324415        Numpar* numpar=NULL;
    46334416
    4634         ParameterInputs* inputs=NULL;
    4635 
    4636         /*recover pointers: */
    4637         inputs=(ParameterInputs*)vinputs;
    4638 
    46394417        /*recover objects from hooks: */
    46404418        nodes=(Node**)hnodes.deliverp();
    … …  
    46434421        numpar=(Numpar*)hnumpar.delivers();
    46444422
     4423        /*retrieve inputs :*/
     4424        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
     4425
    46454426        /*Get out if shelf*/
    4646         if(this->properties.shelf) return;
     4427        if(shelf)return;
    46474428
    46484429        /* Get node coordinates and dof list: */
    … …  
    46504431        GetDofList1(&doflist1[0]);
    46514432
    4652         /* Set grade_g to 0: */
    4653         for(i=0;i<numgrids;i++) grade_g[i]=0.0;
    4654 
    4655         /* recover input parameters: */
    4656         inputs->Recover("drag",&this->properties.k[0],1,dofs1,numgrids,(void**)nodes);
    4657         inputs->Recover("bed",&this->properties.b[0],1,dofs1,numgrids,(void**)nodes);
    4658         inputs->Recover("thickness",&this->properties.h[0],1,dofs1,numgrids,(void**)nodes);
    4659         if(!inputs->Recover("velocity",&vxvy_list[0][0],2,dofs2,numgrids,(void**)nodes)){
    4660                 ISSMERROR("missing velocity input parameter");
    4661         }
    4662         if(!inputs->Recover("adjoint",&adjxadjy_list[0][0],2,dofs2,numgrids,(void**)nodes)){
    4663                 ISSMERROR("missing adjoint input parameter");
    4664         }
    4665 
    4666         /*Initialize parameter lists: */
    4667         for(i=0;i<numgrids;i++){
    4668                 vx_list[i]=vxvy_list[i][0];
    4669                 vy_list[i]=vxvy_list[i][1];
    4670                 adjx_list[i]=adjxadjy_list[i][0];
    4671                 adjy_list[i]=adjxadjy_list[i][1];
    4672         }
     4433        /*Recover inputs: */
     4434        inputs->GetParameterValues(&vx_list[0],&gaussgrids[0][0],3,VxAverageEnum);
     4435        inputs->GetParameterValues(&vy_list[0],&gaussgrids[0][0],3,VyAverageEnum);
     4436        inputs->GetParameterValues(&thickness_list[0],&gaussgrids[0][0],3,ThicknessEnum);
     4437        inputs->GetParameterValues(&bed_list[0],&gaussgrids[0][0],3,BedEnum);
     4438        inputs->GetParameterValues(&dragcoefficient_list[0],&gaussgrids[0][0],3,DragCoefficientEnum);
     4439        inputs->GetParameterValue(&drag_p,DragPEnum);
     4440        inputs->GetParameterValue(&drag_q,DragQEnum);
     4441        inputs->GetParameterValue(&drag_type,DragTypeEnum);
    46734442
    46744443        /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
    … …  
    46844453
    46854454                /*Build alpha_complement_list: */
    4686                 if (!this->properties.shelf && (this->properties.friction_type==2)){
     4455                if (drag_type==2){
    46874456               
    46884457                        /*Allocate friction object: */
    … …  
    46934462                        strcpy(friction->element_type,"2d");
    46944463                       
     4464                       
    46954465                        friction->gravity=matpar->GetG();
    46964466                        friction->rho_ice=matpar->GetRhoIce();
    46974467                        friction->rho_water=matpar->GetRhoWater();
    4698                         friction->K=&this->properties.k[0];
    4699                         friction->bed=&this->properties.b[0];
    4700                         friction->thickness=&this->properties.h[0];
    4701                         friction->velocities=&vxvy_list[0][0];
    4702                         friction->p=this->properties.p;
    4703                         friction->q=this->properties.q;
    4704 
     4468                        friction->K=&dragcoefficient_list[0];
     4469                        friction->bed=&bed_list[0];
     4470                        friction->thickness=&thickness_list[0];
     4471                        friction->vx=&vx_list[0];
     4472                        friction->vy=&vy_list[0];
     4473                        friction->p=drag_p;
     4474                        friction->q=drag_q;
     4475
     4476                       
    47054477                        if(friction->p!=1) ISSMERROR("non-linear friction not supported yet in control methods!");
    47064478
    … …  
    47194491                /*Recover alpha_complement and k: */
    47204492                GetParameterValue(&alpha_complement, &alpha_complement_list[0],gauss_l1l2l3);
    4721                 GetParameterValue(&drag, &this->properties.k[0],gauss_l1l2l3);
     4493                inputs->GetParameterValue(&drag, gauss_l1l2l3,DragCoefficientEnum);
    47224494
    47234495                /*recover lambda and mu: */
    4724                 GetParameterValue(&lambda, &adjx_list[0],gauss_l1l2l3);
    4725                 GetParameterValue(&mu, &adjy_list[0],gauss_l1l2l3);
     4496                inputs->GetParameterValue(&lambda, gauss_l1l2l3,AdjointxEnum);
     4497                inputs->GetParameterValue(&mu, gauss_l1l2l3,AdjointyEnum);
    47264498                       
    47274499                /*recover vx and vy: */
    4728                 GetParameterValue(&vx, &vx_list[0],gauss_l1l2l3);
    4729                 GetParameterValue(&vy, &vy_list[0],gauss_l1l2l3);
     4500                inputs->GetParameterValue(&vx, gauss_l1l2l3,VxEnum);
     4501                inputs->GetParameterValue(&vy, gauss_l1l2l3,VyEnum);
    47304502
    47314503                /* Get Jacobian determinant: */
    … …  
    47394511
    47404512                /*Get k derivative: dk/dx */
    4741                 GetParameterDerivativeValue(&dk[0], &this->properties.k[0],&xyz_list[0][0], gauss_l1l2l3);
     4513                inputs->GetParameterDerivativeValue(&dk[0],&xyz_list[0][0],&gauss_l1l2l3[0],DragCoefficientEnum);
    47424514
    47434515                /*Build gradje_g_gaussian vector (actually -dJ/ddrag): */
    … …  
    47774549/*}}}*/
    47784550/*FUNCTION Tria::GradjDragStokes {{{1*/
    4779 void  Tria::GradjDragStokes(Vec grad_g,void* vinputs,int analysis_type,int sub_analysis_type){
     4551void  Tria::GradjDragStokes(Vec grad_g,int analysis_type,int sub_analysis_type){
    47804552
    47814553        int i;
    … …  
    48374609        double epsilon[3]; /* epsilon=[exx,eyy,exy];*/
    48384610
    4839         ParameterInputs* inputs=NULL;
     4611        /*inputs: */
     4612        bool shelf;
     4613        int  drag_type;
    48404614
    48414615        /*dynamic objects pointed to by hooks: */
    … …  
    48454619        Numpar* numpar=NULL;
    48464620
     4621        /*retrieve inputs :*/
     4622        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
     4623        inputs->GetParameterValue(&drag_type,DragTypeEnum);
     4624
    48474625        /*Get out if shelf*/
    4848         if(this->properties.shelf) return;
    4849 
    4850         /*recover pointers: */
    4851         inputs=(ParameterInputs*)vinputs;
     4626        if(shelf)return;
    48524627
    48534628        /*recover objects from hooks: */
    … …  
    48974672
    48984673                /*Build alpha_complement_list: */
    4899                 if (!this->properties.shelf && (this->properties.friction_type==2)){
     4674                if (drag_type==2){
    49004675
    49014676                        /*Allocate friction object: */
    … …  
    50804855/*}}}*/
    50814856/*FUNCTION Tria::Misfit {{{1*/
    5082 double Tria::Misfit(void* vinputs,int analysis_type,int sub_analysis_type){
     4857double Tria::Misfit(int analysis_type,int sub_analysis_type){
    50834858
    50844859        int i;
    … …  
    51274902        double fit=-1;
    51284903
    5129         ParameterInputs* inputs=NULL;
    5130 
    51314904        /*dynamic objects pointed to by hooks: */
    51324905        Node**  nodes=NULL;
    … …  
    51354908        Numpar* numpar=NULL;
    51364909
     4910        /*inputs: */
     4911        bool onwater;
     4912
     4913        /*retrieve inputs :*/
     4914        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     4915
    51374916        /*If on water, return 0: */
    5138         if(this->properties.onwater)return 0;
    5139 
    5140         /*recover pointers: */
    5141         inputs=(ParameterInputs*)vinputs;
     4917        if(onwater)return 0;
    51424918
    51434919        /*recover objects from hooks: */
    … …  
    53295105/*}}}*/
    53305106/*FUNCTION Tria::SurfaceArea {{{1*/
    5331 double Tria::SurfaceArea(void* vinputs,int analysis_type,int sub_analysis_type){
     5107double Tria::SurfaceArea(int analysis_type,int sub_analysis_type){
    53325108
    53335109        int i;
    … …  
    53465122        Node**  nodes=NULL;
    53475123
     5124        /*inputs: */
     5125        bool onwater;
     5126
     5127        /*retrieve inputs :*/
     5128        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     5129
    53485130        /*recover objects from hooks: */
    53495131        nodes=(Node**)hnodes.deliverp();
    53505132
    53515133        /*If on water, return 0: */
    5352         if(this->properties.onwater)return 0;
     5134        if(onwater)return 0;
    53535135
    53545136        /* Get node coordinates and dof list: */

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

    r3599 r3612  
    99class Element;
    1010class Hook;
    11 class ElementProperties;
    1211class DataSet;
     12class Inputs;
    1313class Node;
    1414struct IoModel;
    … …  
    1818#include "./Hook.h"
    1919#include "./Node.h"
    20 #include "./ElementProperties.h"
    2120#include "../ModelProcessorx/IoModel.h"
    2221#include "../DataSet/DataSet.h"
     22#include "../DataSet/Inputs.h"
    2323
    2424class Tria: public Element{
    … …  
    3333                Hook hnumpar; //hook to 1 numpar
    3434
    35                 ElementProperties properties;
    36                 DataSet* inputs;
     35                Inputs* inputs;
    3736       
    3837        public:
    … …  
    4039                /*FUNCTION constructors, destructors {{{1*/
    4140                Tria();
    42                 Tria(int tria_id,int* tria_node_ids, int tria_matice_id, int tria_matpar_id, int tria_numpar_id, ElementProperties* tria_properties);
    43                 Tria(int tria_id,Hook* tria_hnodes, Hook* tria_hmatice, Hook* tria_hmatpar, Hook* tria_hnumpar, ElementProperties* tria_properties,DataSet* tria_inputs);
     41                Tria(int tria_id,int* tria_node_ids, int tria_matice_id, int tria_matpar_id, int tria_numpar_id);
     42                Tria(int tria_id,Hook* tria_hnodes, Hook* tria_hmatice, Hook* tria_hmatpar, Hook* tria_hnumpar, Inputs* tria_inputs);
    4443                Tria(int i, IoModel* iomodel);
    4544                ~Tria();
    … …  
    6059                /*}}}*/
    6160                /*FUNCTION element numerical routines {{{1*/
    62                 void  CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    63                 void  CreatePVector(Vec pg, void* inputs, int analysis_type,int sub_analysis_type);
     61                void  CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type);
     62                void  CreatePVector(Vec pg, int analysis_type,int sub_analysis_type);
    6463                                void  GetDofList(int* doflist,int* pnumberofdofs);
    6564                void  GetDofList1(int* doflist);
    66                 void  CreateKMatrixDiagnosticHoriz(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    67                 void  CreateKMatrixDiagnosticHorizFriction(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    68                 void  CreateKMatrixDiagnosticSurfaceVert(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    69                 void  CreateKMatrixSlopeCompute(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type);
     65                void  CreateKMatrixDiagnosticHoriz(Mat Kgg,int analysis_type,int sub_analysis_type);
     66                void  CreateKMatrixDiagnosticHorizFriction(Mat Kgg,int analysis_type,int sub_analysis_type);
     67                void  CreateKMatrixDiagnosticSurfaceVert(Mat Kgg,int analysis_type,int sub_analysis_type);
     68                void  CreateKMatrixSlopeCompute(Mat Kgg,int analysis_type,int sub_analysis_type);
    7069                void  GetParameterValue(double* pp, double* plist, double* gauss_l1l2l3);
    7170                void  GetParameterDerivativeValue(double* p, double* plist,double* xyz_list, double* gauss_l1l2l3);
    … …  
    8382                void  GetJacobianInvert(double*  Jinv, double* xyz_list,double* gauss_l1l2l3);
    8483                void  GetJacobian(double* J, double* xyz_list,double* gauss_l1l2l3);
    85                 void  Du(Vec du_g,void* inputs,int analysis_type,int sub_analysis_type);
    86                 void  Gradj(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type,char* control_type);
    87                 void  GradjDrag(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type);
    88                 void  GradjDragStokes(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type);
     84                void  Du(Vec du_g,int analysis_type,int sub_analysis_type);
     85                void  Gradj(Vec grad_g,int analysis_type,int sub_analysis_type,char* control_type);
     86                void  GradjDrag(Vec grad_g,int analysis_type,int sub_analysis_type);
     87                void  GradjDragStokes(Vec grad_g,int analysis_type,int sub_analysis_type);
    8988                void  SurfaceNormal(double* surface_normal, double xyz_list[3][3]);
    90                 void  GradjB(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type);
    91                 double Misfit(void* inputs,int analysis_type,int sub_analysis_type);
    92                 double SurfaceArea(void* inputs,int analysis_type,int sub_analysis_type);
    93                 double CostFunction(void* inputs,int analysis_type,int sub_analysis_type);
     89                void  GradjB(Vec grad_g,int analysis_type,int sub_analysis_type);
     90                double Misfit(int analysis_type,int sub_analysis_type);
     91                double SurfaceArea(int analysis_type,int sub_analysis_type);
     92                double CostFunction(int analysis_type,int sub_analysis_type);
    9493
    95                 void  CreatePVectorDiagnosticHoriz(Vec pg,void* inputs,int analysis_type,int sub_analysis_type);
    96                 void  CreatePVectorDiagnosticBaseVert(Vec pg,void* inputs,int analysis_type,int sub_analysis_type);
    97                 void  CreatePVectorSlopeCompute( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type);
     94                void  CreatePVectorDiagnosticHoriz(Vec pg,int analysis_type,int sub_analysis_type);
     95                void  CreatePVectorDiagnosticBaseVert(Vec pg,int analysis_type,int sub_analysis_type);
     96                void  CreatePVectorSlopeCompute( Vec pg, int analysis_type,int sub_analysis_type);
    9897                void* GetMatPar();
    99                 int   GetShelf();
     98                bool  GetShelf();
    10099                void  GetNodes(void** nodes);
    101                 int   GetOnBed();
    102 
     100                bool  GetOnBed();
    103101                void  GetThicknessList(double* thickness_list);
    104102                void  GetBedList(double* bed_list);
    105                 void  ComputeBasalStress(Vec sigma_b,void* inputs,int analysis_type,int sub_analysis_type);
    106                 void  ComputePressure(Vec p_g,void* inputs,int analysis_type,int sub_analysis_type);
    107                 void  ComputeStrainRate(Vec eps,void* inputs,int analysis_type,int sub_analysis_type);
    108                 void  CreateKMatrixThermal(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    109                 void  CreateKMatrixMelting(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    110                 void  CreatePVectorThermalShelf( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type);
    111                 void  CreatePVectorThermalSheet( Vec pg, void* vinputs, int analysis_type,int sub_analysis_type);
    112                 void  CreateKMatrixPrognostic(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type);
    113                 void  CreatePVectorPrognostic(Vec pg,void* vinputs,int analysis_type,int sub_analysis_type);
    114                 void  CreateKMatrixPrognostic2(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type);
    115                 void  CreatePVectorPrognostic2(Vec pg,void* vinputs,int analysis_type,int sub_analysis_type);
    116                 void  CreateKMatrixBalancedthickness(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type);
    117                 void  CreatePVectorBalancedthickness(Vec pg,void* vinputs,int analysis_type,int sub_analysis_type);
    118                 void  CreateKMatrixBalancedthickness2(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type);
    119                 void  CreatePVectorBalancedthickness2(Vec pg,void* vinputs,int analysis_type,int sub_analysis_type);
    120                 void  CreateKMatrixBalancedvelocities(Mat Kgg,void* vinputs,int analysis_type,int sub_analysis_type);
    121                 void  CreatePVectorBalancedvelocities(Vec pg,void* vinputs,int analysis_type,int sub_analysis_type);
     103                void  ComputeBasalStress(Vec sigma_b,int analysis_type,int sub_analysis_type);
     104                void  ComputePressure(Vec p_g,int analysis_type,int sub_analysis_type);
     105                void  ComputeStrainRate(Vec eps,int analysis_type,int sub_analysis_type);
     106                void  CreateKMatrixThermal(Mat Kgg,int analysis_type,int sub_analysis_type);
     107                void  CreateKMatrixMelting(Mat Kgg,int analysis_type,int sub_analysis_type);
     108                void  CreatePVectorThermalShelf( Vec pg, int analysis_type,int sub_analysis_type);
     109                void  CreatePVectorThermalSheet( Vec pg, int analysis_type,int sub_analysis_type);
     110                void  CreateKMatrixPrognostic(Mat Kgg,int analysis_type,int sub_analysis_type);
     111                void  CreatePVectorPrognostic(Vec pg,int analysis_type,int sub_analysis_type);
     112                void  CreateKMatrixPrognostic2(Mat Kgg,int analysis_type,int sub_analysis_type);
     113                void  CreatePVectorPrognostic2(Vec pg,int analysis_type,int sub_analysis_type);
     114                void  CreateKMatrixBalancedthickness(Mat Kgg,int analysis_type,int sub_analysis_type);
     115                void  CreatePVectorBalancedthickness(Vec pg,int analysis_type,int sub_analysis_type);
     116                void  CreateKMatrixBalancedthickness2(Mat Kgg,int analysis_type,int sub_analysis_type);
     117                void  CreatePVectorBalancedthickness2(Vec pg,int analysis_type,int sub_analysis_type);
     118                void  CreateKMatrixBalancedvelocities(Mat Kgg,int analysis_type,int sub_analysis_type);
     119                void  CreatePVectorBalancedvelocities(Vec pg,int analysis_type,int sub_analysis_type);
    122120                double MassFlux(double* segment,double* ug);
    123121                double GetArea(void);
    … …  
    126124                /*updates:*/
    127125                void  UpdateFromDakota(void* inputs);
    128                 void  UpdateFromInputs(void* inputs);
    129126                void  UpdateInputs(double* solution, int analysis_type, int sub_analysis_type);
    130127                void  UpdateInputsDiagnosticHoriz( double* solution,int analysis_type,int sub_analysis_type);

  • issm/trunk/src/c/objects/TriaVertexInput.h

    r3599 r3612  
    33 */
    44
    5 #include "./Einput.h"
     5#include "./Input.h"
    66
    77#ifndef _TRIAVERTEXINPUT_H_
    88#define _TRIAVERTEXINPUT_H_
    99
    10 class TriaVertexInput: public Einput{
     10class TriaVertexInput: public Input{
    1111
    1212        private:

  • issm/trunk/src/c/objects/Vertex.cpp

    r3567 r3612  
    182182}
    183183/*}}}*/
    184 /*FUNCTION UpdateFromInputs {{{2*/
    185 void  Vertex::UpdateFromInputs(void* vinputs){
    186        
    187         ParameterInputs* inputs=NULL;
    188         Vertex* vertex=NULL;
    189         int dof[1]={0};
    190 
    191         vertex=this;
    192 
    193         /*Recover parameter inputs: */
    194         inputs=(ParameterInputs*)vinputs;
    195 
    196         /*Update internal data if inputs holds new values: */
    197         inputs->Recover("x",&x,1,dof,1,(void**)&vertex);
    198         inputs->Recover("y",&y,1,dof,1,(void**)&vertex);
    199         inputs->Recover("z",&z,1,dof,1,(void**)&vertex);
    200        
    201 }
    202 /*}}}*/
    203184/*FUNCTION UpdateVertexPosition {{{2*/
    204185void  Vertex::UpdatePosition(double* thickness,double* bed){

  • issm/trunk/src/c/objects/Vertex.h

    r3464 r3612  
    4646                int   MyRank();
    4747                void  UpdateFromDakota(void* vinputs);
    48                 void  UpdateFromInputs(void* vinputs);
    4948                void  UpdatePosition(double* thickness,double* bed);
    5049

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

    r3420 r3612  
    2525class NodeSets;
    2626class Model;
     27class TriaVertexInput;
     28class DoubleInput;
     29class IntInput;
     30class BoolInput;
     31class Input;
    2732
    2833/*Abstract class: */
    … …  
    5055#include "./NodeSets.h"
    5156#include "./Model.h"
     57#include "./Input.h"
     58#include "./TriaVertexInput.h"
     59#include "./BoolInput.h"
     60#include "./IntInput.h"
     61#include "./DoubleInput.h"
    5262
    5363/*C objects: */
    5464#include "./Contour.h"
    55 #include "./ParameterInputs.h"
    56 #include "./Input.h"
    5765#include "./Friction.h"
    5866#include "./SolverEnum.h"

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

    r3588 r3612  
    1010
    1111struct OptArgs;
    12 class ParameterInputs;
    1312class FemModel;
    1413
    15 void gradjcompute_core(DataSet* results,Model* model, ParameterInputs* inputs);
     14void gradjcompute_core(DataSet* results,Model* model);
    1615
    17 void diagnostic_core(DataSet* results,Model* model, ParameterInputs* inputs);
    18 void prognostic_core(DataSet* results,Model* model, ParameterInputs* inputs);
    19 void prognostic2_core(DataSet* results,Model* model, ParameterInputs* inputs);
    20 void balancedthickness_core(DataSet* results,Model* model, ParameterInputs* inputs);
    21 void balancedthickness2_core(DataSet* results,Model* model, ParameterInputs* inputs);
    22 void balancedvelocities_core(DataSet* results,Model* model, ParameterInputs* inputs);
    23 void slopecompute_core(DataSet* results,Model* model, ParameterInputs* inputs);
    24 void control_core(DataSet* results,Model* model, ParameterInputs* inputs);
     16void diagnostic_core(DataSet* results,Model* model);
     17void prognostic_core(DataSet* results,Model* model);
     18void prognostic2_core(DataSet* results,Model* model);
     19void balancedthickness_core(DataSet* results,Model* model);
     20void balancedthickness2_core(DataSet* results,Model* model);
     21void balancedvelocities_core(DataSet* results,Model* model);
     22void slopecompute_core(DataSet* results,Model* model);
     23void control_core(DataSet* results,Model* model);
    2524
    26 void thermal_core(DataSet* results,Model* model, ParameterInputs* inputs);
    27 void thermal_core_nonlinear(Vec* ptg,double* pmelting_offset,FemModel* fem,ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
     25void thermal_core(DataSet* results,Model* model);
     26void thermal_core_nonlinear(Vec* ptg,double* pmelting_offset,FemModel* fem,int analysis_type,int sub_analysis_type);
    2827
    29 void steadystate_core(DataSet* results,Model* model, ParameterInputs* inputs);
     28void steadystate_core(DataSet* results,Model* model);
    3029
    31 void diagnostic_core_nonlinear(Vec* pug,Mat* pK_ff0,Mat* pK_fs0, DataSet* loads, FemModel* fem,ParameterInputs* inputs,int analysis_type,int sub_analysis_type);
    32 void diagnostic_core_linear(Vec* ppg,FemModel* fem,ParameterInputs* inputs,int  analysis_type,int sub_analysis_type);
     30void diagnostic_core_nonlinear(Vec* pug,Mat* pK_ff0,Mat* pK_fs0, DataSet* loads, FemModel* fem,int analysis_type,int sub_analysis_type);
     31void diagnostic_core_linear(Vec* ppg,FemModel* fem,int  analysis_type,int sub_analysis_type);
    3332void convergence(int* pconverged, Mat K_ff,Vec p_f,Vec u_f,Vec u_f_old,DataSet* parameters);
    3433
    35 void transient_core(DataSet* results,Model* model, ParameterInputs* inputs);
    36 void transient_core_2d(DataSet* results,Model* model, ParameterInputs* inputs);
    37 void transient_core_3d(DataSet* results,Model* model, ParameterInputs* inputs);
     34void transient_core(DataSet* results,Model* model);
     35void transient_core_2d(DataSet* results,Model* model);
     36void transient_core_3d(DataSet* results,Model* model);
    3837
    3938//int GradJOrth(WorkspaceParams* workspaceparams);
    4039
    41 int GoldenSearch(double* psearch_scalar,double* pJ,double xa, double xb, double tolerance, int maxiter, double fit,double optscal,double (*f)(double*,double,double,FemModel*,ParameterInputs*),FemModel* femmodel);
     40int GoldenSearch(double* psearch_scalar,double* pJ,double xa, double xb, double tolerance, int maxiter, double fit,double optscal,double (*f)(double*,double,double,FemModel*),FemModel* femmodel);
    4241
    43 int BrentSearch(double* psearch_scalar,double* pJ,double xa, double xb, double tolerance, int maxiter, double fit,double optscal,double (*f)(double*,double,double,FemModel*,ParameterInputs*),FemModel* femmodel);
     42int BrentSearch(double* psearch_scalar,double* pJ,double xa, double xb, double tolerance, int maxiter, double fit,double optscal,double (*f)(double*,double,double,FemModel*),FemModel* femmodel);
    4443       
    4544double objectivefunctionC(double search_scalar,OptArgs* optargs);
    … …  
    5251void WriteLockFile(char* filename);
    5352
    54 void ControlInitialization(Model* model, ParameterInputs* inputs);
     53void ControlInitialization(Model* model);
    5554void ControlRestart(Model* model,double* param_g);
    5655

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

    r3582 r3612  
    124124        md.g=0;
    125125        md.yts=0;
    126         md.drag=NaN;
    127126        md.drag_type=0;
    128         md.p=NaN;
    129         md.q=NaN;
    130         md.B=NaN;
    131         md.n=NaN;
     127        md.drag_coefficient=NaN;
     128        md.drag_p=NaN;
     129        md.drag_q=NaN;
     130        md.rheology_B=NaN;
     131        md.rheology_n=NaN;
    132132
    133133        %Geometrical parameters
    … …  
    160160        md.vel_bal=NaN;
    161161        md.vel_obs_raw=NaN;
    162         md.accumulation=NaN;
     162        md.accumulation_rate=NaN;
    163163        md.dhdt=NaN;
    164164        md.geothermalflux=NaN;
    … …  
    237237        md.vel=NaN;
    238238        md.temperature=NaN; %temperature solution vector
    239         md.melting=NaN;
     239        md.melting_rate=NaN;
    240240        md.pressure=NaN;
    241241

  • issm/trunk/src/m/classes/public/marshall.m

    r3582 r3612  
    6464WriteData(fid,md.pressure,'Mat','pressure');
    6565WriteData(fid,md.temperature,'Mat','temperature');
    66 WriteData(fid,md.melting,'Mat','melting');
     66WriteData(fid,md.melting_rate,'Mat','melting_rate');
    6767
    6868
    6969WriteData(fid,md.drag_type,'Integer','drag_type');
    70 WriteData(fid,md.drag,'Mat','drag');
    71 WriteData(fid,md.p,'Mat','p');
    72 WriteData(fid,md.q,'Mat','q');
     70WriteData(fid,md.drag_coefficient,'Mat','drag_coefficient');
     71WriteData(fid,md.drag_p,'Mat','drag_p');
     72WriteData(fid,md.drag_q,'Mat','drag_q');
    7373
    7474WriteData(fid,md.elementoniceshelf,'Mat','elementoniceshelf');
    … …  
    8686WriteData(fid,md.geothermalflux,'Mat','geothermalflux');
    8787WriteData(fid,md.melting,'Mat','melting');
    88 WriteData(fid,md.accumulation,'Mat','accumulation');
     88WriteData(fid,md.accumulation_rate,'Mat','accumulation_rate');
    8989WriteData(fid,md.dhdt,'Mat','dhdt');
    9090
    … …  
    9393WriteData(fid,md.rho_ice,'Scalar','rho_ice');
    9494WriteData(fid,md.g,'Scalar','g');
    95 WriteData(fid,md.B,'Mat','B');
    96 WriteData(fid,md.n,'Mat','n');
     95WriteData(fid,md.rheology_B,'Mat','rheology_B');
     96WriteData(fid,md.rheology_n,'Mat','rheology_n');
    9797
    9898%Control methods

  • issm/trunk/src/m/enum/AirEnum.m

    r3599 r3612  
    77%      macro=AirEnum()
    88
    9 macro=77;
     9macro=80;

  • issm/trunk/src/m/enum/DofVecEnum.m

    r3599 r3612  
    77%      macro=DofVecEnum()
    88
    9 macro=71;
     9macro=74;

  • issm/trunk/src/m/enum/GeographyEnum.m

    r3599 r3612  
    77%      macro=GeographyEnum()
    88
    9 macro=72;
     9macro=75;

  • issm/trunk/src/m/enum/IceEnum.m

    r3599 r3612  
    77%      macro=IceEnum()
    88
    9 macro=76;
     9macro=79;

  • issm/trunk/src/m/enum/IceSheetEnum.m

    r3599 r3612  
    77%      macro=IceSheetEnum()
    88
    9 macro=73;
     9macro=76;

  • issm/trunk/src/m/enum/IceShelfEnum.m

    r3599 r3612  
    77%      macro=IceShelfEnum()
    88
    9 macro=74;
     9macro=77;

  • issm/trunk/src/m/enum/MelangeEnum.m

    r3599 r3612  
    77%      macro=MelangeEnum()
    88
    9 macro=78;
     9macro=81;

  • issm/trunk/src/m/enum/ParamEnum.m

    r3599 r3612  
    77%      macro=ParamEnum()
    88
    9 macro=67;
     9macro=70;

  • issm/trunk/src/m/enum/ResultEnum.m

    r3599 r3612  
    77%      macro=ResultEnum()
    88
    9 macro=68;
     9macro=71;

  • issm/trunk/src/m/enum/RgbEnum.m

    r3599 r3612  
    77%      macro=RgbEnum()
    88
    9 macro=69;
     9macro=72;

  • issm/trunk/src/m/enum/SpcEnum.m

    r3599 r3612  
    77%      macro=SpcEnum()
    88
    9 macro=70;
     9macro=73;

  • issm/trunk/src/m/enum/VxEnum.m

    r3599 r3612  
    77%      macro=VxEnum()
    88
    9 macro=79;
     9macro=82;

  • issm/trunk/src/m/enum/VyEnum.m

    r3599 r3612  
    77%      macro=VyEnum()
    88
    9 macro=80;
     9macro=83;

  • issm/trunk/src/m/enum/WaterEnum.m

    r3599 r3612  
    77%      macro=WaterEnum()
    88
    9 macro=75;
     9macro=78;

  • issm/trunk/src/m/solutions/jpl/diagnostic_core_nonlinear.m

    r3599 r3612  
    1414        converged=0; count=1;
    1515        soln(count).u_g=get(inputs,'velocity',[1 1 (m.parameters.numberofdofspernode>=3) (m.parameters.numberofdofspernode>=3)]); %only keep vz if running with more than 3 dofs per node
     16
    1617        soln(count).u_f= Reducevectorgtof(soln(count).u_g,m.nodesets);
    1718               
    … …  
    5657                %Update elements with new solution
    5758                [m.elements]=UpdateInputs(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,soln(count).u_g,analysis_type,sub_analysis_type);
    58                 error;
    5959               
    6060                %Deal with penalty loads
    … …  
    6363                %penalty constraints
    6464                [loads,constraints_converged,num_unstable_constraints] =PenaltyConstraints( m.elements,m.nodes,m.vertices,loads, m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
     65                error;
    6566               
    6667                displaystring(m.parameters.verbose,'%s%i','      number of unstable constraints: ',num_unstable_constraints);

  • issm/trunk/src/mex/PenaltyConstraints/PenaltyConstraints.cpp

    r3484 r3612  
    5050        PenaltyConstraintsx(&converged, &num_unstable_constraints, elements,nodes,vertices, loads,materials,parameters,inputs,analysis_type,sub_analysis_type);
    5151
     52        elements->Echo();
     53        mexErrMsgTxt(" ");
     54
    5255        /*write output datasets: */
    5356        WriteData(LOADS,loads);
Note: See TracChangeset for help on using the changeset viewer.