Changeset 3621


Ignore:
Timestamp:
04/27/10 07:33:58 (15 years ago)
Author:
Eric.Larour
Message:

Almost there

Location:
issm/trunk/src/c
Files:
2 added
20 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk/src/c/DataSet/DataSet.cpp

    r3612 r3621  
    392392}
    393393/*}}}*/
    394 /*FUNCTION DataSet::FindParam(double* pscalar, char* name){{{1*/
    395 int   DataSet::FindParam(double* pscalar, char* name){
    396        
    397         /*Go through a dataset, and find a Param* object
    398          *which parameter name is "name" : */
    399        
    400         vector<Object*>::iterator object;
    401         Param* param=NULL;
    402 
    403         int found=0;
    404 
    405         for ( object=objects.begin() ; object < objects.end(); object++ ){
    406 
    407                 /*Find param type objects: */
    408                 if((*object)->Enum()==ParamEnum){
    409 
    410                         /*Ok, this object is a parameter, recover it and ask which name it has: */
    411                         param=(Param*)(*object);
    412 
    413                         if (strcmp(param->GetParameterName(),name)==0){
    414                                 /*Ok, this is the one! Recover the value of this parameter: */
    415                                 param->GetParameterValue(pscalar);
    416                                 found=1;
    417                                 break;
    418                         }
    419                 }
    420         }
    421         return found;
    422 }
    423 /*}}}*/
    424 /*FUNCTION DataSet::FindParam(int* pinteger,char* name){{{1*/
    425 int   DataSet::FindParam(int* pinteger,char* name){
    426        
    427        
    428         /*Go through a dataset, and find a Param* object
    429          *which parameter name is "name" : */
    430        
    431         vector<Object*>::iterator object;
    432         Param* param=NULL;
    433 
    434         int found=0;
    435 
    436         for ( object=objects.begin() ; object < objects.end(); object++ ){
    437 
    438                 /*Find param type objects: */
    439                 if((*object)->Enum()==ParamEnum){
    440 
    441                         /*Ok, this object is a parameter, recover it and ask which name it has: */
    442                         param=(Param*)(*object);
    443 
    444                         if (strcmp(param->GetParameterName(),name)==0){
    445                                 /*Ok, this is the one! Recover the value of this parameter: */
    446                                 param->GetParameterValue(pinteger);
    447                                 found=1;
    448                                 break;
    449                         }
    450                 }
    451         }
    452         return found;
    453 }
    454 /*}}}*/
    455 /*FUNCTION DataSet::FindParam(char** pstring,char* name){{{1*/
    456 int   DataSet::FindParam(char** pstring,char* name){
    457        
    458         /*Go through a dataset, and find a Param* object
    459          *which parameter name is "name" : */
    460        
    461         vector<Object*>::iterator object;
    462         Param* param=NULL;
    463 
    464         int found=0;
    465 
    466         for ( object=objects.begin() ; object < objects.end(); object++ ){
    467 
    468                 /*Find param type objects: */
    469                 if((*object)->Enum()==ParamEnum){
    470 
    471                         /*Ok, this object is a parameter, recover it and ask which name it has: */
    472                         param=(Param*)(*object);
    473 
    474                         if (strcmp(param->GetParameterName(),name)==0){
    475                                 /*Ok, this is the one! Recover the value of this parameter: */
    476                                 param->GetParameterValue(pstring);
    477                                 found=1;
    478                                 break;
    479                         }
    480                 }
    481         }
    482         return found;
    483 
    484 }
    485 /*}}}*/
    486 /*FUNCTION DataSet::FindParam(char*** pstringarray,int* pM,char* name){{{1*/
    487 int   DataSet::FindParam(char*** pstringarray,int* pM,char* name){
    488        
    489         /*Go through a dataset, and find a Param* object
    490          *which parameter name is "name" : */
    491        
    492         vector<Object*>::iterator object;
    493         Param* param=NULL;
    494 
    495         int found=0;
    496 
    497         for ( object=objects.begin() ; object < objects.end(); object++ ){
    498 
    499                 /*Find param type objects: */
    500                 if((*object)->Enum()==ParamEnum){
    501 
    502                         /*Ok, this object is a parameter, recover it and ask which name it has: */
    503                         param=(Param*)(*object);
    504 
    505                         if (strcmp(param->GetParameterName(),name)==0){
    506                                 /*Ok, this is the one! Recover the value of this parameter: */
    507                                 param->GetParameterValue(pstringarray);
    508                                 if(pM)*pM=param->GetM();
    509                                 found=1;
    510                                 break;
    511                         }
    512                 }
    513         }
    514         return found;
    515 
    516 }
    517 /*}}}*/
    518 /*FUNCTION DataSet::FindParam(double** pdoublearray,int* pM, int* pN,char* name){{{1*/
    519 int   DataSet::FindParam(double** pdoublearray,int* pM, int* pN,char* name){
    520        
    521         /*Go through a dataset, and find a Param* object
    522          *which parameter name is "name" : */
    523        
    524         vector<Object*>::iterator object;
    525         Param* param=NULL;
    526 
    527         int found=0;
    528 
    529         for ( object=objects.begin() ; object < objects.end(); object++ ){
    530 
    531                 /*Find param type objects: */
    532                 if((*object)->Enum()==ParamEnum){
    533 
    534                         /*Ok, this object is a parameter, recover it and ask which name it has: */
    535                         param=(Param*)(*object);
    536 
    537                         if (strcmp(param->GetParameterName(),name)==0){
    538                                 /*Ok, this is the one! Recover the value of this parameter: */
    539                                 param->GetParameterValue(pdoublearray);
    540                                 if(pM)param->GetM();
    541                                 if(pN)param->GetN();
    542                                 found=1;
    543                                 break;
    544                         }
    545                 }
    546         }
    547         return found;
    548 
    549 }
    550 /*}}}*/
    551 /*FUNCTION DataSet::FindParam(Vec* pvec,char* name){{{1*/
    552 int   DataSet::FindParam(Vec* pvec,char* name){
    553        
    554         /*Go through a dataset, and find a Param* object
    555          *which parameter name is "name" : */
    556        
    557         vector<Object*>::iterator object;
    558         Param* param=NULL;
    559 
    560         int found=0;
    561 
    562         for ( object=objects.begin() ; object < objects.end(); object++ ){
    563 
    564                 /*Find param type objects: */
    565                 if((*object)->Enum()==ParamEnum){
    566 
    567                         /*Ok, this object is a parameter, recover it and ask which name it has: */
    568                         param=(Param*)(*object);
    569 
    570                         if (strcmp(param->GetParameterName(),name)==0){
    571                                 /*Ok, this is the one! Recover the value of this parameter: */
    572                                 param->GetParameterValue(pvec);
    573                                 found=1;
    574                                 break;
    575                         }
    576                 }
    577         }
    578         return found;
    579 
    580 }
    581 /*}}}*/
    582 /*FUNCTION DataSet::FindParamMat* pmat,char* name){{{1*/
    583 int   DataSet::FindParam(Mat* pmat,char* name){
    584        
    585         /*Go through a dataset, and find a Param* object
    586          *which parameter name is "name" : */
    587        
    588         vector<Object*>::iterator object;
    589         Param* param=NULL;
    590 
    591         int found=0;
    592 
    593         for ( object=objects.begin() ; object < objects.end(); object++ ){
    594 
    595                 /*Find param type objects: */
    596                 if((*object)->Enum()==ParamEnum){
    597 
    598                         /*Ok, this object is a parameter, recover it and ask which name it has: */
    599                         param=(Param*)(*object);
    600 
    601                         if (strcmp(param->GetParameterName(),name)==0){
    602                                 /*Ok, this is the one! Recover the value of this parameter: */
    603                                 param->GetParameterValue(pmat);
    604                                 found=1;
    605                                 break;
    606                         }
    607                 }
    608         }
    609         return found;
    610 
    611 }
    612 /*}}}*/
    613 /*FUNCTION DataSet::FindParamObject{{{1*/
    614 Object*   DataSet::FindParamObject(char* name){
    615 
    616         /*Go through a dataset, and find a Param* object
    617          *which parameter name is "name" : */
    618 
    619         vector<Object*>::iterator object;
    620         Param* param=NULL;
    621 
    622         for ( object=objects.begin() ; object < objects.end(); object++ ){
    623 
    624                 /*Find param type objects: */
    625                 if((*object)->Enum()==ParamEnum){
    626 
    627                         /*Ok, this object is a parameter, recover it and ask which name it has: */
    628                         param=(Param*)(*object);
    629 
    630                         if (strcmp(param->GetParameterName(),name)==0){
    631                                 /*Ok, this is the one! Return the object: */
    632                                 return (*object);
    633                         }
    634                 }
    635         }
    636         return NULL;
    637 }
    638 /*}}}*/
    639394/*FUNCTION DataSet::FindResult(Vec* presult,char* name){{{1*/
    640395int   DataSet::FindResult(Vec* presult,char* name){
     
    868623        Load* load=NULL;
    869624        Node* node=NULL;
    870         Numpar* numpar=NULL;
    871625
    872626        for ( object=objects.begin() ; object < objects.end(); object++ ){
     
    886640                        node->Configure(nodes,vertices);
    887641                }
    888                 if((*object)->Enum()==NumparEnum){
    889                         numpar=(Numpar*)(*object);
    890                         numpar->Configure(parameters);
    891                 }
    892 
    893642        }
    894643
  • issm/trunk/src/c/DataSet/DataSet.h

    r3612 r3621  
    4848                int   DeleteObject(int id);
    4949                int   Size();
    50                 int   FindParam(double* pscalar, char* name);
    51                 int   FindParam(int* pinteger,char* name);
    52                 int   FindParam(char** pstring,char* name);
    53                 int   FindParam(char*** pstringarray,int* pM,char* name);
    54                 int   FindParam(double** pdoublearray,int* pM,int* pN,char* name);
    55                 int   FindParam(Vec* pvec,char* name);
    56                 int   FindParam(Mat* pmat,char* name);
     50               
    5751                int   FindResult(Vec* presult,char* name);
    5852                Object* FindParamObject(char* name);
  • issm/trunk/src/c/DataSet/Inputs.cpp

    r3612 r3621  
    249249}
    250250/*}}}*/
     251/*FUNCTION Inputs::GetStrainRateStokes(double* epsilon,double* xyz_list, double* gauss, int xenum, int yenum,int zenum){{{1*/
     252void Inputs::GetStrainRate(double* epsilon,double* xyz_list, double* gauss, int xenum, int yenum,int zenum){
     253
     254        vector<Object*>::iterator object;
     255        Input* xinput=NULL;
     256        Input* yinput=NULL;
     257        Input* yinput=NULL;
     258
     259        /*Go through inputs and find data for xenum: */
     260        for ( object=objects.begin() ; object < objects.end(); object++ ){
     261                xinput=(Input*)(*object);
     262                if (xinput->EnumType()==xenum)break;
     263        }
     264        /*Go through inputs and find data for yenum: */
     265        for ( object=objects.begin() ; object < objects.end(); object++ ){
     266                yinput=(Input*)(*object);
     267                if (yinput->EnumType()==yenum)break;
     268        }
     269        /*Go through inputs and find data for zenum: */
     270        for ( object=objects.begin() ; object < objects.end(); object++ ){
     271                zinput=(Input*)(*object);
     272                if (zinput->EnumType()==zenum)break;
     273        }
     274
     275        if (!xinput | !yinput | !zinput){
     276                /*we could not find one input with the correct enum type. No defaults values were provided,
     277                 * error out: */
     278                ISSMERROR("%s%i%s%i%s%i\n"," could not find input with enum type ",xenum," or enum type ",yenum, " or enum type ",zenum);
     279        }
     280
     281        /*Ok, we have the inputs, call bilinear operator: */
     282        xinput->GetStrainRateStokes(epsilon,yinput,zinput,xyz_list,gauss);
     283
     284}
     285/*}}}*/
    251286/*FUNCTION Inputs::AddInput{{{1*/
    252287int  Inputs::AddInput(Input* in_input){
  • issm/trunk/src/c/DataSet/Inputs.h

    r3612 r3621  
    66#define _INPUTS_H_
    77
     8class Input;
    89#include "./DataSet.h"
    910
     
    2829                void GetParameterDerivativeValue(double* derivativevalues, double* xyz_list, double* gauss,int enum_type);
    2930                void GetStrainRate(double* epsilon,double* xyz_list, double* gauss, int xenum, int yenum);
     31                void GetStrainRateStokes(double* epsilon,double* xyz_list, double* gauss, int xenum, int yenum,int yenum);
    3032                /*}}}*/
    3133
  • issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h

    r3612 r3621  
    9999        InputEnum,
    100100        TriaVertexInputEnum,
     101        SingVertexInputEnum,
     102        BeamVertexInputEnum,
     103        PentaVertexInputEnum,
    101104        BoolInputEnum,
    102105        IntInputEnum,
     
    155158        ElementOnSurfaceEnum,
    156159        SurfaceAreaEnum,
     160        SurfaceSlopexEnum,
     161        SurfaceSlopeyEnum,
    157162        WeightsEnum,
    158163        FitEnum,
    159164        AdjointxEnum,
    160165        AdjointyEnum,
     166        AdjointzEnum,
     167        CollapseEnum,
     168        TemperatureEnum,
    161169        PressureEnum
    162170        /*}}}*/
  • issm/trunk/src/c/Makefile.am

    r3612 r3621  
    6565                                        ./objects/TriaVertexInput.h\
    6666                                        ./objects/TriaVertexInput.cpp\
     67                                        ./objects/SingVertexInput.h\
     68                                        ./objects/SingVertexInput.cpp\
     69                                        ./objects/BeamVertexInput.h\
     70                                        ./objects/BeamVertexInput.cpp\
     71                                        ./objects/PentaVertexInput.h\
     72                                        ./objects/PentaVertexInput.cpp\
    6773                                        ./objects/BoolInput.h\
    6874                                        ./objects/BoolInput.cpp\
     
    8187                                        ./objects/Matpar.h\
    8288                                        ./objects/Matpar.cpp\
    83                                         ./objects/Numpar.h\
    84                                         ./objects/Numpar.cpp\
    8589                                        ./objects/Input.h\
    86                                         ./objects/Input.cpp\
    87                                         ./objects/Einput.h\
    8890                                        ./objects/Spc.cpp\
    8991                                        ./objects/Spc.h\
     
    108110                                        ./DataSet/Inputs.cpp\
    109111                                        ./DataSet/Inputs.h\
     112                                        ./DataSet/Parameters.cpp\
     113                                        ./DataSet/Parameters.h\
    110114                                        ./shared/shared.h\
    111115                                        ./shared/Alloc/alloc.h\
     
    477481                                        ./objects/TriaVertexInput.h\
    478482                                        ./objects/TriaVertexInput.cpp\
     483                                        ./objects/SingVertexInput.h\
     484                                        ./objects/SingVertexInput.cpp\
     485                                        ./objects/BeamVertexInput.h\
     486                                        ./objects/BeamVertexInput.cpp\
     487                                        ./objects/PentaVertexInput.h\
     488                                        ./objects/PentaVertexInput.cpp\
    479489                                        ./objects/BoolInput.h\
    480490                                        ./objects/BoolInput.cpp\
     
    493503                                        ./objects/Matpar.h\
    494504                                        ./objects/Matpar.cpp\
    495                                         ./objects/Numpar.h\
    496                                         ./objects/Numpar.cpp\
    497                                         ./objects/Einput.h\
     505                                        ./objects/Input.h\
    498506                                        ./objects/Spc.cpp\
    499507                                        ./objects/Spc.h\
     
    518526                                        ./DataSet/Inputs.cpp\
    519527                                        ./DataSet/Inputs.h\
     528                                        ./DataSet/Parameters.cpp\
     529                                        ./DataSet/Parameters.h\
    520530                                        ./shared/shared.h\
    521531                                        ./shared/Threads/issm_threads.h\
  • issm/trunk/src/c/ModelProcessorx/CreateDataSets.cpp

    r3588 r3621  
    1515
    1616
    17 void CreateDataSets(DataSet** pelements,DataSet** pnodes, DataSet** pvertices, DataSet** pmaterials, DataSet** pconstraints, DataSet** ploads,DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle){
     17void CreateDataSets(DataSet** pelements,DataSet** pnodes, DataSet** pvertices, DataSet** pmaterials, DataSet** pconstraints, DataSet** ploads,Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle){
    1818
    1919        /*create parameters common to all solutions: */
  • issm/trunk/src/c/ModelProcessorx/CreateParameters.cpp

    r3588 r3621  
    3232        parameters   = new DataSet(ParametersEnum);
    3333
    34         //solution parameters: always 1'st, to speed-up its lookup by elements..
    35         count++;
    36         numpar= new Numpar(count);
    37         parameters->AddObject(numpar);
    38 
    3934        //outputfilename
    4035        count++;
     
    4338        parameters->AddObject(param);
    4439
    45         //analysis and subanalysis
    46         count++;
    47         param= new Param(count,"analysis_type",DOUBLE);
    48         param->SetDouble(iomodel->analysis_type);
    49         parameters->AddObject(param);
    50 
    51         count++;
    52         param= new Param(count,"sub_analysis_type",DOUBLE);
    53         param->SetDouble(iomodel->sub_analysis_type);
    54         parameters->AddObject(param);
    55 
    56        
    5740        //dimension 2d or 3d:
    5841        if (strcmp(iomodel->meshtype,"2d")==0)dim=2;
     
    121104        param->SetDouble(iomodel->dt);
    122105        parameters->AddObject(param);
     106
    123107
    124108        /*ndt: */
  • issm/trunk/src/c/ModelProcessorx/IoModel.h

    r3612 r3621  
    99#include "../io/io.h"
    1010#include "../DataSet/DataSet.h"
     11#include "../DataSet/Parameters.h"
    1112#include "../toolkits/toolkits.h"
    1213
     
    198199
    199200        /*Creation of fem datasets: general drivers*/
    200         void  CreateDataSets(DataSet** pelements,DataSet** pnodes, DataSet** pvertices, DataSet** pmaterials, DataSet** pconstraints, DataSet** ploads,DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    201         void  CreateParameters(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     201        void  CreateDataSets(DataSet** pelements,DataSet** pnodes, DataSet** pvertices, DataSet** pmaterials, DataSet** pconstraints, DataSet** ploads,Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     202        void  CreateParameters(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    202203
    203204       
     
    208209        void    CreateConstraintsDiagnosticHoriz(DataSet** pconstraints,IoModel* iomodel,ConstDataHandle iomodel_handle);
    209210        void  CreateLoadsDiagnosticHoriz(DataSet** ploads, IoModel* iomodel, ConstDataHandle iomodel_handle);
    210         void  CreateParametersDiagnosticHoriz(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     211        void  CreateParametersDiagnosticHoriz(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    211212
    212213        /*diagnostic vertical*/
     
    231232
    232233        /*control:*/
    233         void  CreateParametersControl(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     234        void  CreateParametersControl(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    234235
    235236        /*thermal:*/
     
    237238        void    CreateConstraintsThermal(DataSet** pconstraints,IoModel* iomodel,ConstDataHandle iomodel_handle);
    238239        void  CreateLoadsThermal(DataSet** ploads, IoModel* iomodel, ConstDataHandle iomodel_handle);
    239         void  CreateParametersThermal(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     240        void  CreateParametersThermal(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    240241
    241242        /*melting:*/
     
    243244        void    CreateConstraintsMelting(DataSet** pconstraints,IoModel* iomodel,ConstDataHandle iomodel_handle);
    244245        void  CreateLoadsMelting(DataSet** ploads, IoModel* iomodel, ConstDataHandle iomodel_handle);
    245         void  CreateParametersMelting(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     246        void  CreateParametersMelting(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    246247
    247248        /*prognostic:*/
     
    249250        void    CreateConstraintsPrognostic(DataSet** pconstraints,IoModel* iomodel,ConstDataHandle iomodel_handle);
    250251        void  CreateLoadsPrognostic(DataSet** ploads, IoModel* iomodel, ConstDataHandle iomodel_handle);
    251         void  CreateParametersPrognostic(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     252        void  CreateParametersPrognostic(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    252253
    253254        /*prognostic2:*/
     
    255256        void    CreateConstraintsPrognostic2(DataSet** pconstraints,IoModel* iomodel,ConstDataHandle iomodel_handle);
    256257        void  CreateLoadsPrognostic2(DataSet** ploads, IoModel* iomodel, ConstDataHandle iomodel_handle);
    257         void  CreateParametersPrognostic2(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     258        void  CreateParametersPrognostic2(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    258259
    259260        /*balancedthickness:*/
     
    261262        void    CreateConstraintsBalancedthickness(DataSet** pconstraints,IoModel* iomodel,ConstDataHandle iomodel_handle);
    262263        void  CreateLoadsBalancedthickness(DataSet** ploads, IoModel* iomodel, ConstDataHandle iomodel_handle);
    263         void  CreateParametersBalancedthickness(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     264        void  CreateParametersBalancedthickness(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    264265
    265266        /*balancedthickness2:*/
     
    267268        void    CreateConstraintsBalancedthickness2(DataSet** pconstraints,IoModel* iomodel,ConstDataHandle iomodel_handle);
    268269        void  CreateLoadsBalancedthickness2(DataSet** ploads, IoModel* iomodel, ConstDataHandle iomodel_handle);
    269         void  CreateParametersBalancedthickness2(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     270        void  CreateParametersBalancedthickness2(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    270271
    271272        /*balancedvelocities:*/
     
    273274        void    CreateConstraintsBalancedvelocities(DataSet** pconstraints,IoModel* iomodel,ConstDataHandle iomodel_handle);
    274275        void  CreateLoadsBalancedvelocities(DataSet** ploads, IoModel* iomodel, ConstDataHandle iomodel_handle);
    275         void  CreateParametersBalancedvelocities(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     276        void  CreateParametersBalancedvelocities(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    276277
    277278        /*qmu: */
    278         void CreateParametersQmu(DataSet** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
     279        void CreateParametersQmu(Parameters** pparameters,IoModel* iomodel,ConstDataHandle iomodel_handle);
    279280
    280281       
  • issm/trunk/src/c/io/WriteParams.cpp

    r3570 r3621  
    6363        for(i=0;i<nfields;i++){
    6464
    65                 param=(Param*)parameters->GetObjectByOffset(i+1);
     65                param=(Param*)parameters->GetObjectByOffset(i+1); //skip the numpar object
    6666               
    6767                switch(param->GetType()){
  • issm/trunk/src/c/objects/Beam.cpp

    r3620 r3621  
    2323Beam::Beam(){
    2424        this->inputs=NULL;
     25        this->parameters=NULL;
    2526        return;
    2627}
    2728/*}}}*/
    28 /*FUNCTION Beam::Beam(int id, int* node_ids, int matice_id, int matpar_id, int numpar_id, ElementProperties* properties){{{1*/
    29 Beam::Beam(int beam_id,int* beam_node_ids, int beam_matice_id, int beam_matpar_id, int beam_numpar_id):
     29/*FUNCTION Beam::Beam(int id, int* node_ids, int matice_id, int matpar_id){{{1*/
     30Beam::Beam(int beam_id,int* beam_node_ids, int beam_matice_id, int beam_matpar_id):
    3031        hnodes(beam_node_ids,2),
    3132        hmatice(&beam_matice_id,1),
    32         hmatpar(&beam_matpar_id,1),
    33         hnumpar(&beam_numpar_id,1)
     33        hmatpar(&beam_matpar_id,1)
    3434{
    3535
    3636        /*all the initialization has been done by the initializer, just fill in the id: */
    3737        this->id=beam_id;
     38        this->parameters=NULL;
    3839        this->inputs=new Inputs();
    3940}
    4041/*}}}*/
    41 /*FUNCTION Beam::Beam(int id, Hook* hnodes, Hook* hmatice, Hook* hmatpar, Hook* hnumpar, ElementProperties* properties){{{1*/
    42 Beam::Beam(int beam_id,Hook* beam_hnodes, Hook* beam_hmatice, Hook* beam_hmatpar, Hook* beam_hnumpar, Inputs* beam_inputs):
     42/*FUNCTION Beam::Beam(int id, Hook* hnodes, Hook* hmatice, Hook* hmatpar, Parameters* beam_parameters, ElementProperties* properties){{{1*/
     43Beam::Beam(int beam_id,Hook* beam_hnodes, Hook* beam_hmatice, Hook* beam_hmatpar, Parameters* beam_parameters, Inputs* beam_inputs):
    4344        hnodes(beam_hnodes),
    4445        hmatice(beam_hmatice),
    45         hmatpar(beam_hmatpar),
    46         hnumpar(beam_hnumpar)
     46        hmatpar(beam_hmatpar)
    4747{
    4848
     
    5555                this->inputs=new Inputs();
    5656        }
    57         return;
     57        /*point parameters: */
     58        this->parameters=beam_parameters;
    5859}
    5960/*}}}*/
     
    6869        int   beam_matice_id;
    6970        int   beam_matpar_id;
    70         int   beam_numpar_id;
    7171        int   beam_node_ids[2];
    7272        double nodeinputs[2];
     
    7979        beam_matice_id=index+1; //refers to the corresponding material property card
    8080        beam_matpar_id=iomodel->numberofvertices2d*(iomodel->numlayers-1)+1;//refers to the corresponding matpar property card
    81         beam_numpar_id=1;
    8281        beam_node_ids[0]=index+1;
    8382        beam_node_ids[1]=(int)iomodel->uppernodes[index]; //grid that lays right on top
     
    8685        this->hmatice.Init(&beam_matice_id,1);
    8786        this->hmatpar.Init(&beam_matpar_id,1);
    88         this->hnumpar.Init(&beam_numpar_id,1);
    8987
    9088        //intialize inputs, and add as many inputs per element as requested:
     
    113111        if (iomodel->gridonbed) this->inputs->AddInput(new BoolInput(ElementOnBedEnum,(IssmBool)iomodel->gridonbed[index]));
    114112
     113        //this->parameters: we still can't point to it, it may not even exist. Configure will handle this.
     114        this->parameters=NULL;
     115
     116
    115117}
    116118/*}}}*/
     
    118120Beam::~Beam(){
    119121        delete inputs;
    120         return;
     122        this->parameters=NULL;
    121123}
    122124/*}}}*/
     
    124126/*Object management*/
    125127/*FUNCTION Beam::Configure{{{1*/
    126 void  Beam::Configure(void* ploadsin,void* pnodesin,void* pmaterialsin,void* pparametersin){
     128void  Beam::Configure(DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
    127129
    128130        int i;
    129131       
    130         DataSet* loadsin=NULL;
    131         DataSet* nodesin=NULL;
    132         DataSet* materialsin=NULL;
    133         DataSet* parametersin=NULL;
    134 
    135         /*Recover pointers :*/
    136         loadsin=(DataSet*)ploadsin;
    137         nodesin=(DataSet*)pnodesin;
    138         materialsin=(DataSet*)pmaterialsin;
    139         parametersin=(DataSet*)pparametersin;
    140 
    141132        /*Take care of hooking up all objects for this element, ie links the objects in the hooks to their respective
    142133         * datasets, using internal ids and offsets hidden in hooks: */
     
    144135        hmatice.configure(materialsin);
    145136        hmatpar.configure(materialsin);
    146         hnumpar.configure(parametersin);
     137
     138        /*point parameters to real dataset: */
     139        this->parameters=parametersin;
    147140
    148141}
     
    151144Object* Beam::copy() {
    152145       
    153         return new Beam(this->id,&this->hnodes,&this->hmatice,&this->hmatpar,&this->hnumpar,this->inputs);
     146        return new Beam(this->id,&this->hnodes,&this->hmatice,&this->hmatpar,this->parameters,this->inputs);
    154147
    155148}
     
    163156        hmatice.DeepEcho();
    164157        hmatpar.DeepEcho();
    165         hnumpar.DeepEcho();
     158        printf("   parameters\n");
     159        parameters->DeepEcho();
    166160        printf("   inputs\n");
    167161        inputs->DeepEcho();
     
    187181        hmatice.Demarshall(&marshalled_dataset);
    188182        hmatpar.Demarshall(&marshalled_dataset);
    189         hnumpar.Demarshall(&marshalled_dataset);
    190183
    191184        /*demarshall inputs: */
    192185        inputs=(Inputs*)DataSetDemarshallRaw(&marshalled_dataset);
     186
     187        /*parameters: may not exist even yet, so let Configure handle it: */
     188        this->parameters=NULL;
    193189
    194190        /*return: */
     
    205201        hmatice.Echo();
    206202        hmatpar.Echo();
    207         hnumpar.Echo();
     203        printf("   parameters\n");
     204        parameters->Echo();
    208205        printf("   inputs\n");
    209206        inputs->Echo();
     
    236233        hmatice.Marshall(&marshalled_dataset);
    237234        hmatpar.Marshall(&marshalled_dataset);
    238         hnumpar.Marshall(&marshalled_dataset);
    239235
    240236        /*Marshall inputs: */
     
    243239        memcpy(marshalled_dataset,marshalled_inputs,marshalled_inputs_size*sizeof(char));
    244240        marshalled_dataset+=marshalled_inputs_size;
     241       
     242        /*parameters: don't do anything about it. parameters are marshalled somewhere else!*/
    245243
    246244        xfree((void**)&marshalled_inputs);
     
    257255                +hmatice.MarshallSize()
    258256                +hmatpar.MarshallSize()
    259                 +hnumpar.MarshallSize()
    260257                +inputs->MarshallSize()
    261258                +sizeof(int); //sizeof(int) for enum type
     
    292289        Matpar* matpar=NULL;
    293290        Matice* matice=NULL;
    294         Numpar* numpar=NULL;
    295291       
    296292        /*Get dof list on which we will plug the pressure values: */
     
    301297        matpar=(Matpar*)hmatpar.delivers();
    302298        matice=(Matice*)hmatice.delivers();
    303         numpar=(Numpar*)hnumpar.delivers();
    304299
    305300        /*Get node data: */
     
    454449        Matpar* matpar=NULL;
    455450        Matice* matice=NULL;
    456         Numpar* numpar=NULL;
    457451
    458452        /*recover objects from hooks: */
     
    460454        matpar=(Matpar*)hmatpar.delivers();
    461455        matice=(Matice*)hmatice.delivers();
    462         numpar=(Numpar*)hnumpar.delivers();
    463456
    464457        /*recover doflist: */
  • issm/trunk/src/c/objects/Beam.h

    r3620 r3621  
    1313#include "../ModelProcessorx/IoModel.h"
    1414#include "../DataSet/Inputs.h"
     15#include "../DataSet/DataSet.h"
     16#include "../DataSet/Parameters.h"
    1517#include "./Hook.h"
    1618
     
    3032                Hook hmatice; //hook to 1 matice
    3133                Hook hmatpar; //hook to 1 matpar
    32                 Hook hnumpar; //hook to 1 numpar
     34               
     35                Parameters* parameters; //pointer to solution parameters
    3336                Inputs* inputs;
    3437       
     
    3740                /*constructors, destructors: {{{1*/
    3841                Beam();
    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);
     42                Beam(int beam_id,int* beam_node_ids, int beam_matice_id, int beam_matpar_id);
     43                Beam(int beam_id,Hook* beam_hnodes, Hook* beam_hmatice, Hook* beam_hmatpar, Parameters* beam_parameters, Inputs* beam_inputs);
    4144                Beam(int i, IoModel* iomodel);
    4245                ~Beam();
  • issm/trunk/src/c/objects/Pengrid.cpp

    r3570 r3621  
    1818#include "../include/typedefs.h"
    1919#include "../include/macros.h"
    20                
     20#include "../DataSet/DataSet.h"
     21#include "../DataSet/Inputs.h"
     22       
    2123/*Object constructors and destructor*/
    2224/*FUNCTION Pengrid::constructor {{{1*/
    2325Pengrid::Pengrid(){
    24         return;
    25 }
    26 /*}}}1*/
    27 /*FUNCTION Pengrid::creation {{{1*/
     26        this->inputs=NULL;
     27        this->parameters=NULL;
     28       
     29        /*not active, not zigzagging: */
     30        active=0;
     31        zigzag_counter=0;
     32
     33}
     34/*}}}1*/
     35/*FUNCTION Pengrid::Pengrid(int id, int node_ids int matpar_id){{{1*/
     36Pengrid::Pengrid(int pengrid_id,int pengrid_node_id, int pengrid_matpar_id):
     37        hnode(pengrid_node_ids,1),
     38        hmatice(&pengrid_matice_id,1),
     39        hmatpar(&pengrid_matpar_id,1)
     40{
     41
     42        /*all the initialization has been done by the initializer, just fill in the id: */
     43        this->id=pengrid_id;
     44        this->parameters=NULL;
     45        this->inputs=new Inputs();
     46
     47        /*not active, not zigzagging: */
     48        active=0;
     49        zigzag_counter=0;
     50
     51}
     52/*}}}*/
     53/*FUNCTION Pengrid::Pengrid(int id, Hook* hnodes, Hook* hmatice, Hook* hmatpar, DataSet* parameters, Inputs* pengrid_inputs) {{{1*/
     54Pengrid::Pengrid(int pengrid_id,Hook* pengrid_hnode, Hook* pengrid_hmatpar, Parameters* pengrid_parameters, Inputs* pengrid_inputs):
     55        hnode(pengrid_hnode),
     56        hmatpar(pengrid_hmatpar)
     57{
     58
     59        /*all the initialization has been done by the initializer, just fill in the id: */
     60        this->id=pengrid_id;
     61        if(pengrid_inputs){
     62                this->inputs=(Inputs*)pengrid_inputs->Copy();
     63        }
     64        else{
     65                this->inputs=new Inputs();
     66        }
     67        /*point parameters: */
     68        this->parameters=pengrid_parameters;
     69       
     70        /*not active, not zigzagging: */
     71        active=0;
     72        zigzag_counter=0;
     73
     74}
     75/*}}}*/
    2876Pengrid::Pengrid(int    pengrid_id, int pengrid_node_id,int pengrid_mparid, int pengrid_dof, int pengrid_active, double pengrid_penalty_offset,int pengrid_thermal_steadystate,int pengrid_stabilize_constraints){
    2977       
     
    4290        matpar_offset=UNDEF;
    4391
    44         zigzag_counter=0;
    4592
    4693        return;
    4794}
    4895/*}}}1*/
     96/*FUNCTION Pengrid::Pengrid(int i, IoModel* iomodel){{{1*/
     97Pengrid::Pengrid(int index, IoModel* iomodel){ //i is the element index
     98
     99        int i,j;
     100        int tria_node_ids[3];
     101        int tria_matice_id;
     102        int tria_matpar_id;
     103        double nodeinputs[3];
     104
     105        /*id: */
     106        this->id=index+1;
     107       
     108        /*hooks: */
     109        //go recover node ids, needed to initialize the node hook.
     110        if (iomodel->analysis_type==Prognostic2AnalysisEnum || iomodel->analysis_type==Balancedthickness2AnalysisEnum){
     111                /*Discontinuous Galerkin*/
     112                tria_node_ids[0]=3*index+1;
     113                tria_node_ids[1]=3*index+2;
     114                tria_node_ids[2]=3*index+3;
     115        }
     116        else{
     117                /*Continuous Galerkin*/
     118                for(i=0;i<3;i++){
     119                        tria_node_ids[i]=(int)*(iomodel->elements+3*index+i); //ids for vertices are in the elements array from Matlab
     120                }
     121        }
     122        tria_matice_id=index+1; //refers to the corresponding ice material object
     123        tria_matpar_id=iomodel->numberofelements+1; //refers to the constant material parameters object
     124
     125        this->hnodes.Init(tria_node_ids,3);
     126        this->hmatice.Init(&tria_matice_id,1);
     127        this->hmatpar.Init(&tria_matpar_id,1);
     128
     129        //intialize inputs, and add as many inputs per element as requested:
     130        this->inputs=new Inputs();
     131       
     132        if (iomodel->thickness) {
     133                for(i=0;i<3;i++)nodeinputs[i]=iomodel->thickness[tria_node_ids[i]-1];
     134                this->inputs->AddInput(new PengridVertexInput(ThicknessEnum,nodeinputs));
     135        }
     136        if (iomodel->surface) {
     137                for(i=0;i<3;i++)nodeinputs[i]=iomodel->surface[tria_node_ids[i]-1];
     138                this->inputs->AddInput(new PengridVertexInput(SurfaceEnum,nodeinputs));
     139        }
     140        if (iomodel->bed) {
     141                for(i=0;i<3;i++)nodeinputs[i]=iomodel->bed[tria_node_ids[i]-1];
     142                this->inputs->AddInput(new PengridVertexInput(BedEnum,nodeinputs));
     143        }
     144        if (iomodel->drag_coefficient) {
     145                for(i=0;i<3;i++)nodeinputs[i]=iomodel->drag_coefficient[tria_node_ids[i]-1];
     146                this->inputs->AddInput(new PengridVertexInput(DragCoefficientEnum,nodeinputs));
     147
     148                if (iomodel->drag_p) this->inputs->AddInput(new DoubleInput(DragPEnum,iomodel->drag_p[index]));
     149                if (iomodel->drag_q) this->inputs->AddInput(new DoubleInput(DragQEnum,iomodel->drag_q[index]));
     150                this->inputs->AddInput(new IntInput(DragTypeEnum,iomodel->drag_type));
     151
     152        }
     153        if (iomodel->melting_rate) {
     154                for(i=0;i<3;i++)nodeinputs[i]=iomodel->melting_rate[tria_node_ids[i]-1];
     155                this->inputs->AddInput(new PengridVertexInput(MeltingRateEnum,nodeinputs));
     156        }
     157        if (iomodel->accumulation_rate) {
     158                for(i=0;i<3;i++)nodeinputs[i]=iomodel->accumulation_rate[tria_node_ids[i]-1];
     159                this->inputs->AddInput(new PengridVertexInput(AccumulationRateEnum,nodeinputs));
     160        }
     161        if (iomodel->geothermalflux) {
     162                for(i=0;i<3;i++)nodeinputs[i]=iomodel->geothermalflux[tria_node_ids[i]-1];
     163                this->inputs->AddInput(new PengridVertexInput(GeothermalFluxEnum,nodeinputs));
     164        }       
     165
     166        if (iomodel->elementoniceshelf) this->inputs->AddInput(new BoolInput(ElementOnIceShelfEnum,(IssmBool)iomodel->elementoniceshelf[index]));
     167        if (iomodel->elementonbed) this->inputs->AddInput(new BoolInput(ElementOnBedEnum,(IssmBool)iomodel->elementonbed[index]));
     168        if (iomodel->elementonwater) this->inputs->AddInput(new BoolInput(ElementOnWaterEnum,(IssmBool)iomodel->elementonwater[index]));
     169        if (iomodel->elementonsurface) this->inputs->AddInput(new BoolInput(ElementOnSurfaceEnum,(IssmBool)iomodel->elementonsurface[index]));
     170
     171        //this->parameters: we still can't point to it, it may not even exist. Configure will handle this.
     172        this->parameters=NULL;
     173
     174
     175}
     176/*}}}*/
    49177/*FUNCTION Pengrid::destructor {{{1*/
    50178Pengrid::~Pengrid(){
     
    162290/*FUNCTION Pengrid::CreateKMatrix {{{1*/
    163291
    164 void  Pengrid::CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     292void  Pengrid::CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type){
    165293
    166294        /*No loads applied, do nothing: */
     
    170298/*}}}1*/
    171299/*FUNCTION Pengrid::CreatePVector {{{1*/
    172 void  Pengrid::CreatePVector(Vec pg, void* inputs, int analysis_type,int sub_analysis_type){
     300void  Pengrid::CreatePVector(Vec pg, int analysis_type,int sub_analysis_type){
    173301
    174302        /*No loads applied, do nothing: */
     
    253381/*}}}1*/
    254382/*FUNCTION Pengrid::PenaltyConstrain {{{1*/
    255 void  Pengrid::PenaltyConstrain(int* punstable,void* vinputs,int analysis_type,int sub_analysis_type){
     383void  Pengrid::PenaltyConstrain(int* punstable,int analysis_type,int sub_analysis_type){
    256384
    257385        if ((analysis_type==DiagnosticAnalysisEnum) && ((sub_analysis_type==StokesAnalysisEnum))){
     
    263391        else if (analysis_type==ThermalAnalysisEnum){
    264392               
    265                 PenaltyConstrainThermal(punstable,vinputs,analysis_type,sub_analysis_type);
     393                PenaltyConstrainThermal(punstable,analysis_type,sub_analysis_type);
    266394               
    267395        }
     
    279407/*}}}1*/
    280408/*FUNCTION Pengrid::PenaltyConstrainThermal {{{1*/
    281 void  Pengrid::PenaltyConstrainThermal(int* punstable,void* vinputs,int analysis_type,int sub_analysis_type){
     409void  Pengrid::PenaltyConstrainThermal(int* punstable,int analysis_type,int sub_analysis_type){
    282410
    283411        //   The penalty is stable if it doesn't change during to successive iterations.   
     
    360488/*}}}1*/
    361489/*FUNCTION Pengrid::PenaltyCreateMatrix {{{1*/
    362 void  Pengrid::PenaltyCreateKMatrix(Mat Kgg,void* inputs,double kmax,int analysis_type,int sub_analysis_type){
     490void  Pengrid::PenaltyCreateKMatrix(Mat Kgg,double kmax,int analysis_type,int sub_analysis_type){
    363491
    364492        if ((analysis_type==DiagnosticAnalysisEnum) && ((sub_analysis_type==StokesAnalysisEnum))){
    365493
    366                 PenaltyCreateKMatrixDiagnosticStokes( Kgg,inputs,kmax,analysis_type,sub_analysis_type);
     494                PenaltyCreateKMatrixDiagnosticStokes( Kgg,kmax,analysis_type,sub_analysis_type);
    367495        }
    368496        else if (analysis_type==ThermalAnalysisEnum){
    369497               
    370                 PenaltyCreateKMatrixThermal( Kgg,inputs,kmax,analysis_type,sub_analysis_type);
     498                PenaltyCreateKMatrixThermal( Kgg,kmax,analysis_type,sub_analysis_type);
    371499               
    372500        }
    373501        else if (analysis_type==MeltingAnalysisEnum){
    374502                       
    375                 PenaltyCreateKMatrixMelting( Kgg,inputs,kmax,analysis_type,sub_analysis_type);
     503                PenaltyCreateKMatrixMelting( Kgg,kmax,analysis_type,sub_analysis_type);
    376504
    377505        }
     
    383511/*}}}1*/
    384512/*FUNCTION Pengrid::PenaltyCreateKMatrixDiagnosticStokes {{{1*/
    385 void  Pengrid::PenaltyCreateKMatrixDiagnosticStokes(Mat Kgg,void* vinputs,double kmax,int analysis_type,int sub_analysis_type){
     513void  Pengrid::PenaltyCreateKMatrixDiagnosticStokes(Mat Kgg,double kmax,int analysis_type,int sub_analysis_type){
    386514       
    387515        const int numgrids=1;
     
    421549/*}}}1*/
    422550/*FUNCTION Pengrid::PenaltyCreateKMatrixMelting {{{1*/
    423 void  Pengrid::PenaltyCreateKMatrixMelting(Mat Kgg,void* vinputs,double kmax,int analysis_type,int sub_analysis_type){
     551void  Pengrid::PenaltyCreateKMatrixMelting(Mat Kgg,double kmax,int analysis_type,int sub_analysis_type){
    424552
    425553
     
    469597/*}}}1*/
    470598/*FUNCTION Pengrid::PenaltyCreateKMatrixThermal {{{1*/
    471 void  Pengrid::PenaltyCreateKMatrixThermal(Mat Kgg,void* vinputs,double kmax,int analysis_type,int sub_analysis_type){
     599void  Pengrid::PenaltyCreateKMatrixThermal(Mat Kgg,double kmax,int analysis_type,int sub_analysis_type){
    472600
    473601        int found=0;
     
    498626/*}}}1*/
    499627/*FUNCTION Pengrid::PenaltyCreatePVector {{{1*/
    500 void  Pengrid::PenaltyCreatePVector(Vec pg,void* inputs,double kmax,int analysis_type,int sub_analysis_type){
     628void  Pengrid::PenaltyCreatePVector(Vec pg,double kmax,int analysis_type,int sub_analysis_type){
    501629
    502630        if (analysis_type==ThermalAnalysisEnum){
    503631               
    504                 PenaltyCreatePVectorThermal( pg,inputs,kmax,analysis_type,sub_analysis_type);
     632                PenaltyCreatePVectorThermal( pg,kmax,analysis_type,sub_analysis_type);
    505633               
    506634        }
    507635        else if (analysis_type==MeltingAnalysisEnum){
    508636                       
    509                 PenaltyCreatePVectorMelting( pg,inputs,kmax,analysis_type,sub_analysis_type);
     637                PenaltyCreatePVectorMelting( pg,kmax,analysis_type,sub_analysis_type);
    510638
    511639        }
     
    523651/*}}}1*/
    524652/*FUNCTION Pengrid::PenaltyCreatePVectorMelting {{{1*/
    525 void  Pengrid::PenaltyCreatePVectorMelting(Vec pg, void* vinputs, double kmax,int analysis_type,int sub_analysis_type){
     653void  Pengrid::PenaltyCreatePVectorMelting(Vec pg, double kmax,int analysis_type,int sub_analysis_type){
    526654       
    527655        const int numgrids=1;
     
    595723/*}}}1*/
    596724/*FUNCTION Pengrid::PenaltyCreatePVectorThermal {{{1*/
    597 void  Pengrid::PenaltyCreatePVectorThermal(Vec pg, void* vinputs, double kmax,int analysis_type,int sub_analysis_type){
     725void  Pengrid::PenaltyCreatePVectorThermal(Vec pg, double kmax,int analysis_type,int sub_analysis_type){
    598726
    599727        const int numgrids=1;
  • issm/trunk/src/c/objects/Pengrid.h

    r3463 r3621  
    1616
    1717                int             id;
    18                 int     dof;
    19                 int     active;
    20                 double  penalty_offset;
    21                 int     thermal_steadystate;
    2218               
    23                 /*nodes: */
    24                 int     node_id;
    25                 Node*   node;
    26                 int     node_offset;
     19                Hook hnode;  //hook to 1 node
     20                Hook hmatpar; //hook to 1 matpar
    2721
    28                 int mparid;
    29                 Matpar* matpar;
    30                 int   matpar_offset;
    31 
    32                 int stabilize_constraints;
     22                Parameters* parameters; //pointer to solution parameters
     23                Inputs*  inputs;
     24       
     25                /*internals: */
     26                int active;
    3327                int zigzag_counter;
    3428
    3529        public:
    3630
     31                /*FUNCTION constructors, destructors {{{1*/
    3732                Pengrid();
    38                 Pengrid(int     id, int node_id,int mparid,int dof, int active, double penalty_offset,int thermal_steadystate,int stabilize_constraints);
     33                Pengrid(int pengrid_id,int pengrid_node_id int pengrid_matpar_id);
     34                Pengrid(int pengrid_id,Hook* pengrid_hnode, Hook* pengrid_hmatpar, Parameters* pengrid_parameters, Inputs* pengrid_inputs);
     35                Pengrid(int i, IoModel* iomodel);
    3936                ~Pengrid();
    40 
     37                /*}}}*/
     38                /*FUNCTION object management {{{1*/
     39                void  Configure(void* elements,void* nodes,void* materials);
     40                Object* copy();
     41                void  DeepEcho();
     42                void  Demarshall(char** pmarshalled_dataset);
    4143                void  Echo();
    42                 void  DeepEcho();
     44                int   Enum();
     45                int   GetId();
     46                char* GetName();
    4347                void  Marshall(char** pmarshalled_dataset);
    4448                int   MarshallSize();
    45                 char* GetName();
    46                 void  Demarshall(char** pmarshalled_dataset);
    47                 int   Enum();
    48                 int   GetId();
    4949                int   MyRank();
     50                /*}}}*/
     51                /*FUNCTION element numerical routines {{{1*/
    5052                void  DistributeNumDofs(int* numdofspernode,int analysis_type,int sub_analysis_type);
    51                 void  Configure(void* elements,void* nodes,void* materials);
    52                 void  CreateKMatrix(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type);
    53                 void  CreatePVector(Vec pg, void* inputs, int analysis_type,int sub_analysis_type);
    54                 void  UpdateFromInputs(void* inputs);
    55                 void  PenaltyCreateKMatrix(Mat Kgg,void* inputs,double kmax,int analysis_type,int sub_analysis_type);
    56                 void  PenaltyCreatePVector(Vec pg,void* inputs,double kmax,int analysis_type,int sub_analysis_type);
    57                 void  PenaltyCreateKMatrixDiagnosticStokes(Mat Kgg,void* inputs,double kmax,int analysis_type,int sub_analysis_type);
     53                void  CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type);
     54                void  CreatePVector(Vec pg, int analysis_type,int sub_analysis_type);
     55                void  PenaltyCreateKMatrix(Mat Kgg,double kmax,int analysis_type,int sub_analysis_type);
     56                void  PenaltyCreatePVector(Vec pg,double kmax,int analysis_type,int sub_analysis_type);
     57                void  PenaltyCreateKMatrixDiagnosticStokes(Mat Kgg,double kmax,int analysis_type,int sub_analysis_type);
    5858                void  GetDofList(int* doflist,int* pnumberofdofspernode);
    59                 Object* copy();
    60                 void  PenaltyCreateKMatrixThermal(Mat Kgg,void* vinputs,double kmax,int analysis_type,int sub_analysis_type);
    61                 void  PenaltyCreateKMatrixMelting(Mat Kgg,void* vinputs,double kmax,int analysis_type,int sub_analysis_type);
     59                void  PenaltyCreateKMatrixThermal(Mat Kgg,double kmax,int analysis_type,int sub_analysis_type);
     60                void  PenaltyCreateKMatrixMelting(Mat Kgg,double kmax,int analysis_type,int sub_analysis_type);
    6261                void  MatparConfiguration(Matpar* matpar,int matpar_offset);
    63                 void  PenaltyCreatePVectorThermal(Vec pg, void* inputs, double kmax,int analysis_type,int sub_analysis_type);
    64                 void  PenaltyCreatePVectorMelting(Vec pg, void* inputs, double kmax,int analysis_type,int sub_analysis_type);
    65                 void  PenaltyConstrain(int* punstable,void* inputs,int analysis_type,int sub_analysis_type);
    66                 void  PenaltyConstrainThermal(int* punstable,void* inputs,int analysis_type,int sub_analysis_type);
     62                void  PenaltyCreatePVectorThermal(Vec pg, double kmax,int analysis_type,int sub_analysis_type);
     63                void  PenaltyCreatePVectorMelting(Vec pg, double kmax,int analysis_type,int sub_analysis_type);
     64                void  PenaltyConstrain(int* punstable,int analysis_type,int sub_analysis_type);
     65                void  PenaltyConstrainThermal(int* punstable,int analysis_type,int sub_analysis_type);
     66               
     67                /*updates:*/
     68                void  UpdateFromDakota(void* inputs);
     69                void  UpdateInputs(double* solution, int analysis_type, int sub_analysis_type);
     70                /*}}}*/
    6771
    6872};
  • issm/trunk/src/c/objects/Penta.cpp

    r3620 r3621  
    1717#include "../include/typedefs.h"
    1818#include "../include/macros.h"
     19#include "../DataSet/DataSet.h"
     20#include "../DataSet/Inputs.h"
    1921
    2022/*Object constructors and destructor*/
     
    2224Penta::Penta(){
    2325        this->inputs=NULL;
    24         return;
     26        this->parameters=NULL;
    2527}
    2628/*}}}*/
    2729/*FUNCTION Penta constructor {{{1*/
    28 Penta::Penta(int penta_id,int* penta_node_ids, int penta_matice_id, int penta_matpar_id, int penta_numpar_id):
     30Penta::Penta(int penta_id,int* penta_node_ids, int penta_matice_id, int penta_matpar_id):
    2931        hnodes(penta_node_ids,6),
    3032        hmatice(&penta_matice_id,1),
    31         hmatpar(&penta_matpar_id,1),
    32         hnumpar(&penta_numpar_id,1)
     33        hmatpar(&penta_matpar_id,1)
    3334{
    3435
    3536        /*all the initialization has been done by the initializer, just fill in the id: */
    3637        this->id=penta_id;
     38        this->parameters=NULL;
    3739        this->inputs=new Inputs();
    3840
    3941}
    4042/*}}}*/
    41 /*FUNCTION Penta other constructor {{{1*/
    42 Penta::Penta(int penta_id,Hook* penta_hnodes, Hook* penta_hmatice, Hook* penta_hmatpar, Hook* penta_hnumpar, Inputs* penta_inputs):
     43/*FUNCTION Penta::Penta(int id, Hook* hnodes, Hook* hmatice, Hook* hmatpar, DataSet* parameters, Inputs* penta_inputs) {{{1*/
     44Penta::Penta(int penta_id,Hook* penta_hnodes, Hook* penta_hmatice, Hook* penta_hmatpar, Parameters* penta_parameters, Inputs* penta_inputs):
    4345        hnodes(penta_hnodes),
    4446        hmatice(penta_hmatice),
    45         hmatpar(penta_hmatpar),
    46         hnumpar(penta_hnumpar)
     47        hmatpar(penta_hmatpar)
    4748{
    4849
     
    5556                this->inputs=new Inputs();
    5657        }
     58        /*point parameters: */
     59        this->parameters=penta_parameters;
    5760}
    5861/*}}}*/
     
    6467        int penta_matice_id;
    6568        int penta_matpar_id;
    66         int penta_numpar_id;
    6769        double nodeinputs[6];
    6870        bool collapse;
     
    7779        penta_matice_id=index+1; //refers to the corresponding ice material object
    7880        penta_matpar_id=iomodel->numberofelements+1; //refers to the constant material parameters object
    79         penta_numpar_id=1; //refers to numerical parameters object
    8081
    8182        this->hnodes.Init(&penta_node_ids[0],6);
    8283        this->hmatice.Init(&penta_matice_id,1);
    8384        this->hmatpar.Init(&penta_matpar_id,1);
    84         this->hnumpar.Init(&penta_numpar_id,1);
    8585
    8686        //intialize inputs, and add as many inputs per element as requested:
     
    142142Penta::~Penta(){
    143143        delete inputs;
    144         return;
     144        this->parameters=NULL;
    145145}
    146146/*}}}*/
    147147
    148148/*Object management: */
    149 /*FUNCTION Configure {{{1*/
    150 void  Penta::Configure(void* ploadsin,void* pnodesin,void* pmaterialsin,void* pparametersin){
     149/*FUNCTION Penta::Configure {{{1*/
     150void  Penta::Configure(DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
    151151
    152152        int i;
    153 
    154         DataSet* loadsin=NULL;
    155         DataSet* nodesin=NULL;
    156         DataSet* materialsin=NULL;
    157         DataSet* parametersin=NULL;
    158 
    159         /*Recover pointers :*/
    160         loadsin=(DataSet*)ploadsin;
    161         nodesin=(DataSet*)pnodesin;
    162         materialsin=(DataSet*)pmaterialsin;
    163         parametersin=(DataSet*)pparametersin;
    164153
    165154        /*Take care of hooking up all objects for this element, ie links the objects in the hooks to their respective
     
    168157        hmatice.configure(materialsin);
    169158        hmatpar.configure(materialsin);
    170         hnumpar.configure(parametersin);
     159
     160        /*point parameters to real dataset: */
     161        this->parameters=parametersin;
    171162
    172163}
     
    174165/*FUNCTION copy {{{1*/
    175166Object* Penta::copy() {
    176         return new Penta(this->id,&this->hnodes,&this->hmatice,&this->hmatpar,&this->hnumpar,this->inputs);
     167        return new Penta(this->id,&this->hnodes,&this->hmatice,&this->hmatpar,this->parameters,this->inputs);
    177168}
    178169/*}}}*/
     
    195186        hmatice.Demarshall(&marshalled_dataset);
    196187        hmatpar.Demarshall(&marshalled_dataset);
    197         hnumpar.Demarshall(&marshalled_dataset);
    198188
    199189        /*demarshall inputs: */
    200190        inputs=(Inputs*)DataSetDemarshallRaw(&marshalled_dataset);
     191
     192        /*parameters: may not exist even yet, so let Configure handle it: */
     193        this->parameters=NULL;
    201194
    202195        /*return: */
     
    214207        hmatice.DeepEcho();
    215208        hmatpar.DeepEcho();
    216         hnumpar.DeepEcho();
     209        printf("   parameters\n");
     210        parameters->DeepEcho();
    217211        printf("   inputs\n");
    218212        inputs->DeepEcho();
     
    230224        hmatice.Echo();
    231225        hmatpar.Echo();
    232         hnumpar.Echo();
     226        printf("   parameters\n");
     227        parameters->Echo();
    233228        printf("   inputs\n");
    234229        inputs->Echo();
     
    267262        hmatice.Marshall(&marshalled_dataset);
    268263        hmatpar.Marshall(&marshalled_dataset);
    269         hnumpar.Marshall(&marshalled_dataset);
    270264
    271265        /*Marshall inputs: */
     
    275269        marshalled_dataset+=marshalled_inputs_size;
    276270
     271        /*parameters: don't do anything about it. parameters are marshalled somewhere else!*/
     272
    277273        xfree((void**)&marshalled_inputs);
    278 
    279274
    280275        *pmarshalled_dataset=marshalled_dataset;
     
    289284                +hmatice.MarshallSize()
    290285                +hmatpar.MarshallSize()
    291                 +hnumpar.MarshallSize()
    292286                +inputs->MarshallSize()
    293287                +sizeof(int); //sizeof(int) for enum type
     
    304298        Hook* tria_hmatice=NULL;
    305299        Hook* tria_hmatpar=NULL;
    306         Hook* tria_hnumpar=NULL;
     300        Parameters* tria_parameters=NULL;
    307301        Inputs* tria_inputs=NULL;
    308302
     
    314308        tria_hmatice=this->hmatice.Spawn(&zero,1);
    315309        tria_hmatpar=this->hmatpar.Spawn(&zero,1);
    316         tria_hnumpar=this->hnumpar.Spawn(&zero,1);
     310        tria_parameters=this->parameters;
    317311        tria_inputs=(Inputs*)this->inputs->Spawn(indices,3);
    318312
    319         tria=new Tria(this->id,tria_hnodes,tria_hmatice,tria_hmatpar,tria_hnumpar,tria_inputs);
     313        tria=new Tria(this->id,tria_hnodes,tria_hmatice,tria_hmatpar,tria_parameters,tria_inputs);
    320314
    321315        delete tria_hnodes;
    322316        delete tria_hmatice;
    323317        delete tria_hmatpar;
    324         delete tria_hnumpar;
    325318        delete tria_inputs;
    326319
     
    497490        bool onbed;
    498491
     492        /*parameters: */
     493        double stokesreconditioning;
     494
    499495        /*dynamic objects pointed to by hooks: */
    500496        Node**  nodes=NULL;
    501497        Matpar* matpar=NULL;
    502498        Matice* matice=NULL;
    503         Numpar* numpar=NULL;
    504499
    505500        /*Check analysis_types*/
     
    510505        matpar=(Matpar*)hmatpar.delivers();
    511506        matice=(Matice*)hmatice.delivers();
    512         numpar=(Numpar*)hnumpar.delivers();
    513507
    514508        /*recover some inputs: */
    515509        inputs->GetParameterValue(&onbed,ElementOnBedEnum);
    516510
     511        /*retrieve some parameters: */
     512        this->parameters->FindParam(&stokesreconditioning,"stokesreconditioning");
     513       
    517514        if(!onbed){
    518515                //put zero
     
    552549
    553550                        /*Compute Stress*/
    554                         sigma_xx=viscosity*epsilon[0]-pressure*numpar->stokesreconditioning; // sigma = nu eps - pressure
    555                         sigma_yy=viscosity*epsilon[1]-pressure*numpar->stokesreconditioning;
    556                         sigma_zz=viscosity*epsilon[2]-pressure*numpar->stokesreconditioning;
     551                        sigma_xx=viscosity*epsilon[0]-pressure*stokesreconditioning; // sigma = nu eps - pressure
     552                        sigma_yy=viscosity*epsilon[1]-pressure*stokesreconditioning;
     553                        sigma_zz=viscosity*epsilon[2]-pressure*stokesreconditioning;
    557554                        sigma_xy=viscosity*epsilon[3];
    558555                        sigma_xz=viscosity*epsilon[4];
     
    723720        else if (analysis_type==BalancedthicknessAnalysisEnum){
    724721
    725                 CreateKMatrixBalancedthickness( Kgg,inputs,analysis_type,sub_analysis_type);
     722                CreateKMatrixBalancedthickness( Kgg,analysis_type,sub_analysis_type);
    726723        }
    727724        else if (analysis_type==BalancedvelocitiesAnalysisEnum){
    728725
    729                 CreateKMatrixBalancedvelocities( Kgg,inputs,analysis_type,sub_analysis_type);
     726                CreateKMatrixBalancedvelocities( Kgg,analysis_type,sub_analysis_type);
    730727        }
    731728        else if (analysis_type==ThermalAnalysisEnum){
     
    745742/*FUNCTION CreateKMatrixBalancedthickness {{{1*/
    746743
    747 void  Penta::CreateKMatrixBalancedthickness(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     744void  Penta::CreateKMatrixBalancedthickness(Mat Kgg,int analysis_type,int sub_analysis_type){
    748745
    749746        /*Collapsed formulation: */
    750747        Tria*  tria=NULL;
    751748
     749        /*flags: */
     750        bool onwater;
     751        bool onbed;
     752
     753        /*recover some inputs: */
     754        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     755        inputs->GetParameterValue(&onbed,ElementOnBedEnum);
     756
    752757        /*If on water, skip: */
    753         if(this->properties.onwater)return;
     758        if(onwater)return;
    754759
    755760        /*Is this element on the bed? :*/
    756         if(!this->properties.onbed)return;
     761        if(!onbed)return;
    757762
    758763        /*Spawn Tria element from the base of the Penta: */
    759764        tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria.
    760         tria->CreateKMatrix(Kgg,inputs, analysis_type,sub_analysis_type);
     765        tria->CreateKMatrix(Kgg, analysis_type,sub_analysis_type);
    761766        delete tria;
    762767        return;
     
    766771/*FUNCTION CreateKMatrixBalancedvelocities {{{1*/
    767772
    768 void  Penta::CreateKMatrixBalancedvelocities(Mat Kgg,void* inputs,int analysis_type,int sub_analysis_type){
     773void  Penta::CreateKMatrixBalancedvelocities(Mat Kgg,int analysis_type,int sub_analysis_type){
    769774
    770775        /*Collapsed formulation: */
    771776        Tria*  tria=NULL;
    772777
     778        /*flags: */
     779        bool onbed;
     780        bool onwater;
     781
     782        /*recover some inputs: */
     783        inputs->GetParameterValue(&onbed,ElementOnBedEnum);
     784        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     785
    773786        /*If on water, skip: */
    774         if(this->properties.onwater)return;
     787        if(onwater)return;
    775788
    776789        /*Is this element on the bed? :*/
    777         if(!this->properties.onbed)return;
     790        if(!onbed)return;
    778791
    779792        /*Spawn Tria element from the base of the Penta: */
    780793        tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria.
    781         tria->CreateKMatrix(Kgg,inputs, analysis_type,sub_analysis_type);
     794        tria->CreateKMatrix(Kgg,analysis_type,sub_analysis_type);
    782795        delete tria;
    783796        return;
     
    860873        double MOUNTAINKEXPONENT=10;
    861874
     875        /*parameters: */
     876        double viscosity_overshoot;
     877
    862878        /*Collapsed formulation: */
    863879        Tria*  tria=NULL;
     
    866882        Node**  nodes=NULL;
    867883        Matice* matice=NULL;
    868         Numpar* numpar=NULL;
    869884
    870885        /*recover objects from hooks: */
    871886        nodes=(Node**)hnodes.deliverp();
    872887        matice=(Matice*)hmatice.delivers();
    873         numpar=(Numpar*)hnumpar.delivers();
    874888
    875889        /*inputs: */
     
    884898        inputs->GetParameterValue(&onbed,ElementOnBedEnum);
    885899        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
     900
     901        /*retrieve some parameters: */
     902        this->parameters->FindParam(&viscosity_overshoot,"viscosity_overshoot");
    886903
    887904        /*If on water, skip stiffness: */
     
    961978                                  onto this scalar matrix, so that we win some computational time: */
    962979
    963                                 newviscosity=viscosity+numpar->viscosity_overshoot*(viscosity-oldviscosity);
     980                                newviscosity=viscosity+viscosity_overshoot*(viscosity-oldviscosity);
    964981                                D_scalar=newviscosity*gauss_weight*Jdet;
    965982                                for (i=0;i<5;i++){
     
    10941111        double  drag_p,drag_q;
    10951112
     1113        /*parameters: */
     1114        double stokesreconditioning;
     1115
    10961116        /*dynamic objects pointed to by hooks: */
    10971117        Node**  nodes=NULL;
    10981118        Matpar* matpar=NULL;
    10991119        Matice* matice=NULL;
    1100         Numpar* numpar=NULL;
    11011120
    11021121        /*inputs: */
     
    11171136        matpar=(Matpar*)hmatpar.delivers();
    11181137        matice=(Matice*)hmatice.delivers();
    1119         numpar=(Numpar*)hnumpar.delivers();
    11201138
    11211139
     
    11241142        rho_ice=matpar->GetRhoIce();
    11251143        gravity=matpar->GetG();
     1144
     1145        /*retrieve some parameters: */
     1146        this->parameters->FindParam(&stokesreconditioning,"stokesreconditioning");
    11261147
    11271148        /* Get node coordinates and dof list: */
     
    11731194                        }
    11741195                        for (i=6;i<8;i++){
    1175                                 D[i][i]=-D_scalar*numpar->stokesreconditioning;
     1196                                D[i][i]=-D_scalar*stokesreconditioning;
    11761197                        }
    11771198
     
    12801301                        DLStokes[9][8]=-viscosity*gauss_weight*Jdet2d*bed_normal[0]/2.0;
    12811302                        DLStokes[10][10]=-viscosity*gauss_weight*Jdet2d*bed_normal[1]/2.0;
    1282                         DLStokes[11][11]=numpar->stokesreconditioning*gauss_weight*Jdet2d*bed_normal[0];
    1283                         DLStokes[12][12]=numpar->stokesreconditioning*gauss_weight*Jdet2d*bed_normal[1];
    1284                         DLStokes[13][13]=numpar->stokesreconditioning*gauss_weight*Jdet2d*bed_normal[2];
     1303                        DLStokes[11][11]=stokesreconditioning*gauss_weight*Jdet2d*bed_normal[0];
     1304                        DLStokes[12][12]=stokesreconditioning*gauss_weight*Jdet2d*bed_normal[1];
     1305                        DLStokes[13][13]=stokesreconditioning*gauss_weight*Jdet2d*bed_normal[2];
    12851306
    12861307                        /*  Do the triple product tL*D*L: */
     
    15701591
    15711592        int     dofs[3]={0,1,2};
    1572         double  dt;
    15731593        double  K[2][2]={0.0};
    15741594
     
    16051625        double     mixed_layer_capacity,thermal_exchange_velocity;
    16061626
     1627        /*parameters: */
     1628        double dt,artdiff,epsvel;
     1629
    16071630        /*dynamic objects pointed to by hooks: */
    16081631        Node**  nodes=NULL;
    16091632        Matpar* matpar=NULL;
    1610         Numpar* numpar=NULL;
    16111633
    16121634        /*Collapsed formulation: */
     
    16301652        nodes=(Node**)hnodes.deliverp();
    16311653        matpar=(Matpar*)hmatpar.delivers();
    1632         numpar=(Numpar*)hnumpar.delivers();
    16331654
    16341655        /* Get node coordinates and dof list: */
     
    16421663        heatcapacity=matpar->GetHeatCapacity();
    16431664        thermalconductivity=matpar->GetThermalConductivity();
    1644         dt=numpar->dt;
     1665
     1666        /*retrieve some parameters: */
     1667        this->parameters->FindParam(&dt,"dt");
     1668        this->parameters->FindParam(&artdiff,"art_diff");
     1669        this->parameters->FindParam(&epsvel,"epsvel");
    16451670
    16461671        /* Get gaussian points and weights. Penta is an extrusion of a Tria, we therefore
     
    17351760
    17361761                        /*Artifficial diffusivity*/
    1737                         if(numpar->artdiff){
     1762                        if(artdiff){
    17381763                                /*Build K: */
    1739                                 D_scalar=gauss_weight*Jdet/(pow(u,2)+pow(v,2)+numpar->epsvel);
     1764                                D_scalar=gauss_weight*Jdet/(pow(u,2)+pow(v,2)+epsvel);
    17401765                                if(dt){
    17411766                                        D_scalar=D_scalar*dt;
     
    18221847        else if (analysis_type==BalancedthicknessAnalysisEnum){
    18231848
    1824                 CreatePVectorPrognostic( pg,inputs,analysis_type,sub_analysis_type);
     1849                CreatePVectorPrognostic( pg,analysis_type,sub_analysis_type);
    18251850        }
    18261851        else if (analysis_type==BalancedvelocitiesAnalysisEnum){
    18271852
    1828                 CreatePVectorPrognostic( pg,inputs,analysis_type,sub_analysis_type);
     1853                CreatePVectorPrognostic( pg,analysis_type,sub_analysis_type);
    18291854        }
    18301855        else if (analysis_type==ThermalAnalysisEnum){
     
    18431868/*}}}*/
    18441869/*FUNCTION CreatePVectorBalancedthickness {{{1*/
    1845 
    1846 void Penta::CreatePVectorBalancedthickness( Vec pg, void* inputs, int analysis_type,int sub_analysis_type){
     1870void Penta::CreatePVectorBalancedthickness( Vec pg, int analysis_type,int sub_analysis_type){
    18471871
    18481872        /*Collapsed formulation: */
    18491873        Tria*  tria=NULL;
    18501874
     1875        /*flags: */
     1876        bool onbed;
     1877        bool onwater;
     1878
     1879        /*recover some inputs: */
     1880        inputs->GetParameterValue(&onbed,ElementOnBedEnum);
     1881        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     1882
    18511883        /*If on water, skip: */
    1852         if(this->properties.onwater)return;
     1884        if(onwater)return;
    18531885
    18541886        /*Is this element on the bed? :*/
    1855         if(!this->properties.onbed)return;
     1887        if(!onbed)return;
    18561888
    18571889        /*Spawn Tria element from the base of the Penta: */
    18581890        tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria.
    1859         tria->CreatePVector(pg,inputs, analysis_type,sub_analysis_type);
     1891        tria->CreatePVector(pg, analysis_type,sub_analysis_type);
    18601892        delete tria;
    18611893        return;
     
    18631895/*}}}*/
    18641896/*FUNCTION CreatePVectorBalancedvelocities {{{1*/
    1865 
    1866 void Penta::CreatePVectorBalancedvelocities( Vec pg, void* inputs, int analysis_type,int sub_analysis_type){
     1897void Penta::CreatePVectorBalancedvelocities( Vec pg, int analysis_type,int sub_analysis_type){
    18671898
    18681899        /*Collapsed formulation: */
    18691900        Tria*  tria=NULL;
    18701901
     1902        /*flags: */
     1903        bool onbed;
     1904        bool onwater;
     1905
     1906        /*recover some inputs: */
     1907        inputs->GetParameterValue(&onbed,ElementOnBedEnum);
     1908        inputs->GetParameterValue(&onwater,ElementOnWaterEnum);
     1909
    18711910        /*If on water, skip: */
    1872         if(this->properties.onwater)return;
     1911        if(onwater)return;
    18731912
    18741913        /*Is this element on the bed? :*/
    1875         if(!this->properties.onbed)return;
     1914        if(!onbed)return;
    18761915
    18771916        /*Spawn Tria element from the base of the Penta: */
    18781917        tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria.
    1879         tria->CreatePVector(pg,inputs, analysis_type,sub_analysis_type);
     1918        tria->CreatePVector(pg, analysis_type,sub_analysis_type);
    18801919        delete tria;
    18811920        return;
     
    21102149        Matpar* matpar=NULL;
    21112150        Matice* matice=NULL;
    2112         Numpar* numpar=NULL;
     2151
     2152        /*parameters: */
     2153        double stokesreconditioning;
    21132154
    21142155        /*inputs: */
     
    21222163        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
    21232164
     2165        /*retrieve some parameters: */
     2166        this->parameters->FindParam(&stokesreconditioning,"stokesreconditioning");
     2167
    21242168        /*If on water, skip load: */
    21252169        if(onwater)return;
     
    21292173        matpar=(Matpar*)hmatpar.delivers();
    21302174        matice=(Matice*)hmatice.delivers();
    2131         numpar=(Numpar*)hnumpar.delivers();
    21322175
    21332176
     
    22002243                        }
    22012244                        for (i=6;i<8;i++){
    2202                                 D[i][i]=-D_scalar*numpar->stokesreconditioning;
     2245                                D[i][i]=-D_scalar*stokesreconditioning;
    22032246                        }
    22042247
     
    25132556        int       num_vert_gauss=3;
    25142557
    2515         double dt;
    25162558        double temperature_list[numgrids];
    25172559        double temperature;
     
    25522594        Tria*  tria=NULL;
    25532595
     2596        /*parameters: */
     2597        double dt;
     2598
    25542599        /*dynamic objects pointed to by hooks: */
    25552600        Node**  nodes=NULL;
    25562601        Matpar* matpar=NULL;
    25572602        Matice* matice=NULL;
    2558         Numpar* numpar=NULL;
    25592603
    25602604        /*inputs: */
     
    25682612        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
    25692613
     2614        /*retrieve some parameters: */
     2615        this->parameters->FindParam(&dt,"dt");
     2616
    25702617        /*If on water, skip: */
    25712618        if(onwater)return;
     
    25752622        matpar=(Matpar*)hmatpar.delivers();
    25762623        matice=(Matice*)hmatice.delivers();
    2577         numpar=(Numpar*)hnumpar.delivers();
    25782624
    25792625        /* Get node coordinates and dof list: */
     
    25882634        beta=matpar->GetBeta();
    25892635        meltingpoint=matpar->GetMeltingPoint();
    2590         dt=numpar->dt;
    25912636
    25922637        /* Get gaussian points and weights. Penta is an extrusion of a Tria, we therefore
  • issm/trunk/src/c/objects/Penta.h

    r3620 r3621  
    2121#include "../ModelProcessorx/IoModel.h"
    2222#include "./Node.h"
     23#include "../DataSet/DataSet.h"
     24#include "../DataSet/Parameters.h"
     25#include "../DataSet/Inputs.h"
    2326
    2427class Penta: public Element{
     
    3033                Hook hmatice; //hook to 1 matice
    3134                Hook hmatpar; //hook to 1 matpar
    32                 Hook hnumpar; //hook to 1 numpar
    3335
     36                Parameters* parameters; //pointer to solution parameters
    3437                Inputs* inputs;
    3538
     
    3841                /*FUNCTION constructors, destructors {{{1*/
    3942                Penta();
    40                 Penta(int penta_id,int* penta_node_ids, int penta_matice_id, int penta_matpar_id, int penta_numpar_id);
    41                 Penta(int penta_id,Hook* penta_hnodes, Hook* penta_hmatice, Hook* penta_hmatpar, Hook* penta_hnumpar, Inputs* inputs);
     43                Penta(int penta_id,int* penta_node_ids, int penta_matice_id, int penta_matpar_id);
     44                Penta(int penta_id,Hook* penta_hnodes, Hook* penta_hmatice, Hook* penta_hmatpar, Parameters* penta_parameters, Inputs* inputs);
    4245                Penta(int i, IoModel* iomodel);
    4346                ~Penta();
     
    104107                void  CreateKMatrixPrognostic(Mat Kgg,int analysis_type,int sub_analysis_type);
    105108                void  CreatePVectorPrognostic( Vec pg,  int analysis_type,int sub_analysis_type);
    106 
     109                void  CreateKMatrixBalancedthickness(Mat Kgg,int analysis_type,int sub_analysis_type);
     110                void  CreateKMatrixBalancedvelocities(Mat Kgg,int analysis_type,int sub_analysis_type);
    107111                void  CreateKMatrixDiagnosticStokes( Mat Kgg,  int analysis_type,int sub_analysis_type);
     112                void  CreatePVectorBalancedthickness( Vec pg, int analysis_type,int sub_analysis_type);
     113                void  CreatePVectorBalancedvelocities( Vec pg, int analysis_type,int sub_analysis_type);
    108114                void  CreatePVectorDiagnosticStokes( Vec pg, int analysis_type,int sub_analysis_type);
    109115                void  ReduceMatrixStokes(double* Ke_reduced, double* Ke_temp);
  • issm/trunk/src/c/objects/Sing.cpp

    r3620 r3621  
    1616#include "../shared/shared.h"
    1717#include "../DataSet/DataSet.h"
     18#include "../DataSet/Parameters.h"
    1819#include "../DataSet/Inputs.h"
    1920#include "../include/typedefs.h"
     
    2425Sing::Sing(){
    2526        this->inputs=NULL;
     27        this->parameters=NULL;
    2628        return;
    2729}
    2830/*}}}*/
    2931/*FUNCTION Sing constructor {{{1*/
    30 Sing::Sing(int sing_id,int* sing_node_ids, int sing_matice_id, int sing_matpar_id, int sing_numpar_id):
     32Sing::Sing(int sing_id,int* sing_node_ids, int sing_matice_id, int sing_matpar_id):
    3133        hnodes(sing_node_ids,1),
    3234        hmatice(&sing_matice_id,1),
    33         hmatpar(&sing_matpar_id,1),
    34         hnumpar(&sing_numpar_id,1)
     35        hmatpar(&sing_matpar_id,1)
    3536{
    3637
    3738        /*all the initialization has been done by the initializer, just fill in the id: */
    3839        this->id=sing_id;
     40        this->parameters=NULL;
    3941        this->inputs=new Inputs();
    4042
     
    4345/*}}}*/
    4446/*FUNCTION Sing other constructor {{{1*/
    45 Sing::Sing(int sing_id,Hook* sing_hnodes, Hook* sing_hmatice, Hook* sing_hmatpar, Hook* sing_hnumpar, Inputs* sing_inputs):
     47Sing::Sing(int sing_id,Hook* sing_hnodes, Hook* sing_hmatice, Hook* sing_hmatpar, Parameters* sing_parameters,Inputs* sing_inputs):
    4648        hnodes(sing_hnodes),
    4749        hmatice(sing_hmatice),
    48         hmatpar(sing_hmatpar),
    49         hnumpar(sing_hnumpar)
     50        hmatpar(sing_hmatpar)
    5051{
    5152
     
    5859                this->inputs=new Inputs();
    5960        }
    60         return;
     61        /*point parameters: */
     62        this->parameters=sing_parameters;
    6163}
    6264/*}}}*/
     
    6668        int sing_matice_id;
    6769        int sing_matpar_id;
    68         int sing_numpar_id;
    6970       
    7071        int sing_g;
     
    7980        sing_matice_id=i+1; //refers to the corresponding material property card
    8081        sing_matpar_id=iomodel->numberofvertices+1;//refers to the corresponding matpar property card
    81         sing_numpar_id=1;
    8282        sing_g=i+1;
    8383
     
    8585        this->hmatice.Init(&sing_matice_id,1);
    8686        this->hmatpar.Init(&sing_matpar_id,1);
    87         this->hnumpar.Init(&sing_numpar_id,1);
    8887
    8988        //intialize inputs, and add as many inputs per element as requested:
     
    9291        if (iomodel->thickness) this->inputs->AddInput(new SingVertexInput(ThicknessEnum,iomodel->thickness[i]));
    9392        if (iomodel->drag_coefficient) this->inputs->AddInput(new SingVertexInput(DragCoefficientEnum,iomodel->drag_coefficient[i]));
     93
     94        //this->parameters: we still can't point to it, it may not even exist. Configure will handle this.
     95        this->parameters=NULL;
    9496
    9597}
     
    98100Sing::~Sing(){
    99101        delete inputs;
    100         return;
     102        this->parameters=NULL;
    101103}
    102104/*}}}*/
     
    104106/*Object management*/
    105107/*FUNCTION Sing::Configure {{{1*/
    106 void  Sing::Configure(void* ploadsin, void* pnodesin,void* pmaterialsin,void* pparametersin){
     108void  Sing::Configure(DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
    107109
    108110        int i;
    109111       
    110         DataSet* loadsin=NULL;
    111         DataSet* nodesin=NULL;
    112         DataSet* materialsin=NULL;
    113         DataSet* parametersin=NULL;
    114 
    115         /*Recover pointers :*/
    116         loadsin=(DataSet*)ploadsin;
    117         nodesin=(DataSet*)pnodesin;
    118         materialsin=(DataSet*)pmaterialsin;
    119         parametersin=(DataSet*)pparametersin;
    120 
    121112        /*Take care of hooking up all objects for this element, ie links the objects in the hooks to their respective
    122113         * datasets, using internal ids and offsets hidden in hooks: */
     
    124115        hmatice.configure(materialsin);
    125116        hmatpar.configure(materialsin);
    126         hnumpar.configure(parametersin);
     117
     118        /*point parameters to real dataset: */
     119        this->parameters=parametersin;
    127120
    128121}
     
    131124Object* Sing::copy() {
    132125       
    133         return new Sing(*this);
     126        return new Sing(this->id,&this->hnodes,&this->hmatice,&this->hmatpar,this->parameters,this->inputs);
    134127
    135128}
     
    143136        hmatice.DeepEcho();
    144137        hmatpar.DeepEcho();
    145         hnumpar.DeepEcho();
     138        printf("   parameters\n");
     139        parameters->DeepEcho();
    146140        printf("   inputs\n");
    147141        inputs->DeepEcho();
     
    167161        hmatice.Demarshall(&marshalled_dataset);
    168162        hmatpar.Demarshall(&marshalled_dataset);
    169         hnumpar.Demarshall(&marshalled_dataset);
    170 
    171         /*demarshall inputs: */
    172         inputs=(Inputs*)DataSetDemarshallRaw(&marshalled_dataset);
     163
     164        /*parameters: may not exist even yet, so let Configure handle it: */
     165        this->parameters=NULL;
    173166
    174167        /*return: */
     
    186179        hmatice.Echo();
    187180        hmatpar.Echo();
    188         hnumpar.Echo();
     181        printf("   parameters\n");
     182        parameters->Echo();
    189183        printf("   inputs\n");
    190184        inputs->Echo();
     
    215209        hmatice.Marshall(&marshalled_dataset);
    216210        hmatpar.Marshall(&marshalled_dataset);
    217         hnumpar.Marshall(&marshalled_dataset);
    218211
    219212        /*Marshall inputs: */
     
    222215        memcpy(marshalled_dataset,marshalled_inputs,marshalled_inputs_size*sizeof(char));
    223216        marshalled_dataset+=marshalled_inputs_size;
     217
     218        /*parameters: don't do anything about it. parameters are marshalled somewhere else!*/
     219
    224220       
    225221        xfree((void**)&marshalled_inputs);
     
    236232                +hmatice.MarshallSize()
    237233                +hmatpar.MarshallSize()
    238                 +hnumpar.MarshallSize()
    239                 +inputs->MarshallSize()
    240234                +sizeof(int); //sizeof(int) for enum type
    241235}
     
    369363        Matpar* matpar=NULL;
    370364        Matice* matice=NULL;
    371         Numpar* numpar=NULL;
    372365
    373366        /*recover objects from hooks: */
     
    375368        matpar=(Matpar*)hmatpar.delivers();
    376369        matice=(Matice*)hmatice.delivers();
    377         numpar=(Numpar*)hnumpar.delivers();
    378370
    379371        inputs->GetParameterValue(&slope[0],SurfaceSlopexEnum);
  • issm/trunk/src/c/objects/Sing.h

    r3620 r3621  
    1313#include "./Hook.h"
    1414#include "../DataSet/Inputs.h"
     15#include "../DataSet/DataSet.h"
     16#include "../DataSet/Parameters.h"
    1517
    1618class Sing: public Element{
     
    2426                Hook hmatice; //hook to 1 matice
    2527                Hook hmatpar; //hook to 1 matpar
    26                 Hook hnumpar; //hook to 1 numpar
     28               
     29                Parameters* parameters; //pointer to solution parameters
    2730                Inputs* inputs;
    2831
     
    3134                /*constructors, destructors: {{{1*/
    3235                Sing();
    33                 Sing(int sing_id,int* sing_node_ids, int sing_matice_id, int sing_matpar_id, int sing_numpar_id);
    34                 Sing(int sing_id,Hook* sing_hnodes, Hook* sing_hmatice, Hook* sing_hmatpar, Hook* sing_hnumpar, Inputs* sing_inputs);
     36                Sing(int sing_id,int* sing_node_ids, int sing_matice_id, int sing_matpar_id);
     37                Sing(int sing_id,Hook* sing_hnodes, Hook* sing_hmatice, Hook* sing_hmatpar, Parameters* sing_parameters,Inputs* sing_inputs);
    3538                Sing(int i, IoModel* iomodel);
    3639                ~Sing();
  • issm/trunk/src/c/objects/Tria.cpp

    r3620 r3621  
    2727Tria::Tria(){
    2828        this->inputs=NULL;
    29         return;
     29        this->parameters=NULL;
    3030}
    3131/*}}}*/
     
    148148        delete inputs;
    149149        this->parameters=NULL;
    150         return;
    151150}
    152151/*}}}*/
  • issm/trunk/src/c/objects/Tria.h

    r3620 r3621  
    66#define _TRIA_H_
    77
    8 class Object;
    9 class Element;
    10 class Hook;
    11 class DataSet;
    12 class Inputs;
    13 class Node;
    14 struct IoModel;
    15 
    168#include "./Object.h"
    179#include "./Element.h"
     
    2012#include "../ModelProcessorx/IoModel.h"
    2113#include "../DataSet/DataSet.h"
     14#include "../DataSet/Parameters.h"
    2215#include "../DataSet/Inputs.h"
    2316
Note: See TracChangeset for help on using the changeset viewer.