Changeset 9356


Ignore:
Timestamp:
08/19/11 18:04:03 (14 years ago)
Author:
Eric.Larour
Message:

Major rewrite of the code so that IoModel now has a parameters dataset, which gets
loaded with all the int,char and double objects in the input file.
This is a lot more flexible, as anyone can add a field to the model, and it will
automatically appear in the IoModel parameters dataset.

Not debugged with respect to nightly runs yet. Trying to get the whole change finished.

Location:
issm/trunk/src/c
Files:
104 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/src/c/Container/Parameters.cpp

    r9320 r9356  
    430430}
    431431/*}}}*/
     432/*FUNCTION Parameters::UnitConversion(int direction_enum);{{{1*/
     433void   Parameters::UnitConversion(int direction_enum){
     434
     435        vector<Object*>::iterator object;
     436        Param* param=NULL;
     437
     438        for ( object=objects.begin() ; object < objects.end(); object++ ){
     439                param=(Param*)(*object);
     440                param->UnitConversion(direction_enum);
     441        }
     442
     443}
     444/*}}}*/
    432445
    433446/*FUNCTION Parameters::FindParamObject{{{1*/

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

    r8600 r9356  
    5454                void  SetParam(Mat mat,int enum_type);
    5555                void  SetParam(FILE* fid,int enum_type);
     56                void  UnitConversion(int direction_enum);
    5657
    5758                Object* FindParamObject(int enum_type);

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

    r9291 r9356  
    177177        AdjointzEnum,
    178178        AdjointpEnum,
    179         ArtDiffEnum,
     179        Fake35Enum,
    180180        BedEnum,
    181181        BathymetryEnum,
    … …  
    187187        Fake2Enum,
    188188        ConstantEnum,
    189         NumControlsEnum,
     189        Fake37Enum,
    190190        ControlTypeEnum,
    191191        ConvergedEnum,
    … …  
    400400        QmuInNameEnum,
    401401        QmuMassFluxSegmentsEnum,
    402         QmuNPartEnum,
     402        Fake39Enum,
    403403        QmuOutNameEnum,
    404404        QmuPartEnum,
    … …  
    422422        FsetEnum,
    423423        SsetEnum,
    424         GroundingLineMigrationEnum,
     424        GroundinglineMigrationEnum,
    425425        YtsEnum,
    426426        /*}}}*/
    … …  
    458458        VelAbsGradientEnum,
    459459        DatasetInputEnum,
    460         NumResponsesEnum,
     460        Fake38Enum,
    461461        StepResponsesEnum,
    462462        IntMatParamEnum,
    … …  
    489489        NodeOnHutterEnum,
    490490        ZEnum,
    491         GlMigrationEnum,
     491        Fake36Enum,
    492492        RiftinfoEnum,
    493493        ElementOnIceSheetEnum,

  • issm/trunk/src/c/EnumDefinitions/EnumToModelField.cpp

    r8926 r9356  
    2222                case VyEnum : return "vy";
    2323                case VyObsEnum : return "vy_obs";
    24                 case GroundingLineMigrationEnum : return "gl_migration";
    2524                case BasalMeltingRateEnum : return "basal_melting_rate";
    2625      case SurfaceAccumulationRateEnum : return "surface_accumulation_rate";

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

    r6413 r9356  
    1111#define SQRT3 1.732050807568877293527446341505872366942805253810380628055806979
    1212#define PI 3.141592653589793238462643383279502884197169399375105820974944592308
     13#define YTS 365.0*24.0*3600.0
    1314
    1415#define NDOF1 1

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

    r9320 r9356  
    2424        int     dummy;
    2525
    26         int     qmu_npart;
     26        int     npart;
    2727        double *qmu_part  = NULL;
    2828        int     numberofvertices;
    … …  
    3939        /*Some parameters: */
    4040        numberofvertices=vertices->NumberOfVertices();
    41         parameters->FindParam(&qmu_npart,QmuNPartEnum);
     41        parameters->FindParam(&npart,NpartEnum);
    4242
    4343        /*average onto the separate areas. The result will be a npart sized vector. */
    4444
    4545        /*allocate: */
    46         partition_contributions=NewVec(qmu_npart);
    47         partition_areas=NewVec(qmu_npart);
    48         vec_average=NewVec(qmu_npart);
     46        partition_contributions=NewVec(npart);
     47        partition_areas=NewVec(npart);
     48        vec_average=NewVec(npart);
    4949
    5050        /*loop on each element, and add contribution of the element to the partition (surface weighted average): */

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

    r8649 r9356  
    2424       
    2525        /*Recover parameters*/
    26         parameters->FindParam(&num_responses,NumResponsesEnum);
     26        parameters->FindParam(&num_responses,NumCmResponsesEnum);
    2727        parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum);
    2828

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

    r9320 r9356  
    3636
    3737        /*retrieve npart: */
    38         parameters->FindParam(&npart,QmuNPartEnum);
     38        parameters->FindParam(&npart,NpartEnum);
    3939
    4040        /*save the d_responses pointer: */

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

    r9291 r9356  
    145145                case AdjointzEnum : return "Adjointz";
    146146                case AdjointpEnum : return "Adjointp";
    147                 case ArtDiffEnum : return "ArtDiff";
     147                case Fake35Enum : return "Fake35";
    148148                case BedEnum : return "Bed";
    149149                case BathymetryEnum : return "Bathymetry";
    … …  
    155155                case Fake2Enum : return "Fake2";
    156156                case ConstantEnum : return "Constant";
    157                 case NumControlsEnum : return "NumControls";
     157                case Fake37Enum : return "Fake37";
    158158                case ControlTypeEnum : return "ControlType";
    159159                case ConvergedEnum : return "Converged";
    … …  
    350350                case QmuInNameEnum : return "QmuInName";
    351351                case QmuMassFluxSegmentsEnum : return "QmuMassFluxSegments";
    352                 case QmuNPartEnum : return "QmuNPart";
     352                case Fake39Enum : return "Fake39";
    353353                case QmuOutNameEnum : return "QmuOutName";
    354354                case QmuPartEnum : return "QmuPart";
    … …  
    372372                case FsetEnum : return "Fset";
    373373                case SsetEnum : return "Sset";
    374                 case GroundingLineMigrationEnum : return "GroundingLineMigration";
     374                case GroundinglineMigrationEnum : return "GroundinglineMigration";
    375375                case YtsEnum : return "Yts";
    376376                case TriangleInterpEnum : return "TriangleInterp";
    … …  
    401401                case VelAbsGradientEnum : return "VelAbsGradient";
    402402                case DatasetInputEnum : return "DatasetInput";
    403                 case NumResponsesEnum : return "NumResponses";
     403                case Fake38Enum : return "Fake38";
    404404                case StepResponsesEnum : return "StepResponses";
    405405                case IntMatParamEnum : return "IntMatParam";
    … …  
    431431                case NodeOnHutterEnum : return "NodeOnHutter";
    432432                case ZEnum : return "Z";
    433                 case GlMigrationEnum : return "GlMigration";
     433                case Fake36Enum : return "Fake36";
    434434                case RiftinfoEnum : return "Riftinfo";
    435435                case ElementOnIceSheetEnum : return "ElementOnIceSheet";

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

    r9211 r9356  
    1919       
    2020        /*retrieve parameters: */
    21         parameters->FindParam(&migration_style,GroundingLineMigrationEnum);
     21        parameters->FindParam(&migration_style,GroundinglineMigrationEnum);
    2222
    2323        /*call different migration modules, according to user wishes: */

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

    r9339 r9356  
    1616       
    1717        int     numberofvertices;
    18         int     qmu_npart;
     18        int     npart;
    1919        double *qmu_part  = NULL;
    2020
    … …  
    2525
    2626        /*retrieve parameters: */
    27         parameters->FindParam(&qmu_npart,QmuNPartEnum);
     27        parameters->FindParam(&npart,NpartEnum);
    2828        parameters->FindParam(&qmu_part,&dummy,QmuPartEnum);
    2929        numberofvertices=vertices->NumberOfVertices();
    … …  
    4040                       
    4141                        /*Variable is scaled. Determine root name of variable (ex: scaled_DragCoefficient_1 -> DragCoefficient). Allocate distributed_values and fill the
    42                          * distributed_values with the next qmu_npart variables: */
     42                         * distributed_values with the next npart variables: */
    4343                       
    4444                        //strcpy(root,strstr(descriptor,"_")+1); *strstr(root,"_")='\0';
    … …  
    4747
    4848
    49                         distributed_values=(double*)xmalloc(qmu_npart*sizeof(double));
    50                         for(j=0;j<qmu_npart;j++){
     49                        distributed_values=(double*)xmalloc(npart*sizeof(double));
     50                        for(j=0;j<npart;j++){
    5151                                distributed_values[j]=variables[i+j];
    5252                        }
    … …  
    7676
    7777                        /*increment i to skip the distributed values just collected: */
    78                         i+=qmu_npart-1; //careful, the for loop will add 1.
     78                        i+=npart-1; //careful, the for loop will add 1.
    7979
    8080                        /*Free allocations: */

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

    r9340 r9356  
    2626        bool    spcpresent=false;
    2727        int     count=0;
     28        int     numberofvertices;
    2829
    2930        /*variables being fetched: */
    … …  
    3132        int     M,N;
    3233
     34        /*Fetch parameters: */
     35        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    3336
    3437        /*First of, find the record for the enum, and get code  of data type: */
    … …  
    4245
    4346        /*Transient or static?:*/
    44         if(M==iomodel->numberofvertices){
     47        if(M==numberofvertices){
    4548                /*static: just create Constraints objects*/
    4649                count=0;
    4750       
    4851                /*Create Constraints from x,y,z: */
    49                 for (i=0;i<iomodel->numberofvertices;i++){
     52                for (i=0;i<numberofvertices;i++){
    5053                       
    5154                        /*keep only this partition's nodes:*/
    … …  
    6063                }
    6164        }
    62         else if (M==iomodel->numberofvertices+1){
     65        else if (M==(numberofvertices+1)){
    6366                /*transient: create transient SpcTransient objects. Same logic, except we need to retrieve
    6467                 * various times and values to initialize an SpcTransient object: */
    … …  
    7477
    7578                /*Create constraints from x,y,z: */
    76                 for (i=0;i<iomodel->numberofvertices;i++){
     79                for (i=0;i<numberofvertices;i++){
    7780                       
    7881                        /*keep only this partition's nodes:*/

  • issm/trunk/src/c/modules/ModelProcessorx/Balancethickness/CreateConstraintsBalancethickness.cpp

    r9340 r9356  
    1010void    CreateConstraintsBalancethickness(Constraints** pconstraints, IoModel* iomodel){
    1111
     12        int    prognostic_DG;   
     13       
     14        /*Fetch parameters: */
     15        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
     16
    1217        /*Output*/
    1318        Constraints* constraints = NULL;
    … …  
    2025
    2126        /*Do not add constraints in DG*/
    22         if(!iomodel->prognostic_DG){
     27        if(!prognostic_DG){
    2328                IoModelToConstraintsx(constraints,iomodel,SpcthicknessEnum,BalancethicknessAnalysisEnum);
    2429        }

  • issm/trunk/src/c/modules/ModelProcessorx/Balancethickness/CreateLoadsBalancethickness.cpp

    r9340 r9356  
    1616        int i;
    1717        int element;
     18        int prognostic_DG;
     19        int numberofedges;
     20
     21        /*Fetch parameters: */
     22        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
    1823
    1924        /*Output*/
    … …  
    2732       
    2833        /*Loads only in DG*/
    29         if (iomodel->prognostic_DG){
     34        if (prognostic_DG){
    3035
    3136                /*Get edges and elements*/
    32                 iomodel->FetchData(&iomodel->edges,&iomodel->numberofedges,NULL,EdgesEnum);
     37                iomodel->FetchData(&iomodel->edges,&numberofedges,NULL,EdgesEnum);
    3338                iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
    3439                iomodel->FetchData(&iomodel->thickness,NULL,NULL,ThicknessEnum);
    3540
    3641                /*First load data:*/
    37                 for (i=0;i<iomodel->numberofedges;i++){
     42                for (i=0;i<numberofedges;i++){
    3843
    3944                        /*Get left and right elements*/

  • issm/trunk/src/c/modules/ModelProcessorx/Balancethickness/CreateNodesBalancethickness.cpp

    r9340 r9356  
    2121        int  io_index;
    2222        bool continuous_galerkin=true;
     23        int    dim;
     24        int    numberofelements;
     25        int    numberofvertices;
     26        int    prognostic_DG;
    2327
    2428        /*DataSets: */
    2529        Nodes*    nodes = NULL;
     30
     31        /*Fetch parameters: */
     32        iomodel->parameters->FindParam(&dim,DimEnum);
     33        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     34        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     35        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
    2636
    2737        /*Recover pointer: */
    … …  
    3242
    3343        /*Continuous Galerkin partition of nodes: */
    34         if(iomodel->prognostic_DG) continuous_galerkin=false;
     44        if(prognostic_DG) continuous_galerkin=false;
    3545        NodesPartitioning(&iomodel->my_nodes,iomodel->my_elements,iomodel->my_vertices,iomodel,continuous_galerkin);
    3646
    3747        /*Check in 3d*/
    38         if(iomodel->prognostic_DG && iomodel->dim==3) _error_("DG 3d not implemented yet");
     48        if(prognostic_DG && dim==3) _error_("DG 3d not implemented yet");
    3949
    4050        /*First fetch data: */
    … …  
    5060
    5161                /*Build Nodes dataset (Continuous Galerkin)*/
    52                 for (i=0;i<iomodel->numberofvertices;i++){
     62                for (i=0;i<numberofvertices;i++){
    5363
    5464                        if(iomodel->my_vertices[i]){
    … …  
    6373
    6474                /*Build Nodes dataset -> 3 for each element (Discontinuous Galerkin)*/
    65                 for (i=0;i<iomodel->numberofelements;i++){
     75                for (i=0;i<numberofelements;i++){
    6676                        for (j=0;j<3;j++){
    6777
    … …  
    7181                                        vertex_id=(int)*(iomodel->elements+3*i+j); //(Matlab indexing)
    7282                                        io_index=vertex_id-1;                      //(C indexing)
    73                                         _assert_(vertex_id>0 && vertex_id<=iomodel->numberofvertices);
     83                                        _assert_(vertex_id>0 && vertex_id<=numberofvertices);
    7484
    7585                                        //Compute Node id

  • issm/trunk/src/c/modules/ModelProcessorx/Balancethickness/UpdateElementsBalancethickness.cpp

    r9343 r9356  
    1616void    UpdateElementsBalancethickness(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    dim;
     19        int    numberofelements;
     20
    1821        /*Fetch data needed: */
     22        iomodel->parameters->FindParam(&dim,DimEnum);
    1923        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
     24        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
    2025
    2126        /*Update elements: */
    2227        int counter=0;
    23         for(int i=0;i<iomodel->numberofelements;i++){
     28        for(int i=0;i<numberofelements;i++){
    2429                if(iomodel->my_elements[i]){
    2530                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    4247        iomodel->FetchDataToInput(elements,DhdtEnum);
    4348
    44         if (iomodel->dim==3){
     49        if (dim==3){
    4550                iomodel->FetchDataToInput(elements,ElementOnBedEnum);
    4651                iomodel->FetchDataToInput(elements,ElementOnSurfaceEnum);

  • issm/trunk/src/c/modules/ModelProcessorx/BedSlope/CreateNodesBedSlope.cpp

    r9340 r9356  
    1818        int i;
    1919        bool continuous_galerkin=true;
     20        int    numberofvertices;
    2021
    2122        /*DataSets: */
    2223        Nodes*    nodes = NULL;
     24
     25        /*Fetch parameters: */
     26        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    2327
    2428        /*Recover pointer: */
    … …  
    3943        iomodel->FetchData(&iomodel->nodeonwater,NULL,NULL,NodeOnWaterEnum);
    4044
    41         for (i=0;i<iomodel->numberofvertices;i++){
     45        for (i=0;i<numberofvertices;i++){
    4246
    4347                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/BedSlope/UpdateElementsBedSlope.cpp

    r9343 r9356  
    1616void    UpdateElementsBedSlope(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    dim;
     19        int    numberofelements;
     20
    1821        /*Fetch data needed: */
     22        iomodel->parameters->FindParam(&dim,DimEnum);
    1923        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
     24        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
    2025
    2126        /*Update elements: */
    2227        int counter=0;
    23         for(int i=0;i<iomodel->numberofelements;i++){
     28        for(int i=0;i<numberofelements;i++){
    2429                if(iomodel->my_elements[i]){
    2530                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    3338        iomodel->FetchDataToInput(elements,ElementOnWaterEnum);
    3439
    35         if (iomodel->dim==3){
     40        if (dim==3){
    3641                iomodel->FetchDataToInput(elements,ElementOnBedEnum);
    3742                iomodel->FetchDataToInput(elements,ElementOnSurfaceEnum);

  • issm/trunk/src/c/modules/ModelProcessorx/Control/CreateParametersControl.cpp

    r9340 r9356  
    1414void CreateParametersControl(Parameters** pparameters,IoModel* iomodel,int solution_type,int analysis_type){
    1515       
    16         int i;
     16        int         i;
    1717        Parameters* parameters=NULL;
     18        bool        control_analysis=false;
     19        int         nsteps;
     20        int         num_control_type;
     21        int         num_cm_responses;
     22        int*        control_type=NULL;
    1823
    1924        /*Get parameters: */
    2025        parameters=*pparameters;
    21         parameters->AddObject(new   BoolParam(ControlAnalysisEnum,iomodel->control_analysis));
    2226
    23         if(iomodel->control_analysis){
     27        /*retrieve some parameters: */
     28        parameters->FindParam(&control_analysis,ControlAnalysisEnum);
    2429
    25                 /*How many controls and how many responses?*/
    26                 parameters->AddObject(new IntParam(NumControlsEnum,iomodel->num_control_type));
    27                 parameters->AddObject(new IntParam(NumResponsesEnum,iomodel->num_cm_responses));
     30        if(control_analysis){
    2831
    2932                /*What control type?*/
    30                 iomodel->FetchData(&iomodel->control_type,NULL,NULL,ControlTypeEnum);
    31                 parameters->AddObject(new IntVecParam(ControlTypeEnum,iomodel->control_type,iomodel->num_control_type));
    32                 xfree((void**)&iomodel->control_type);
     33                iomodel->FetchData(&control_type,&num_control_type,NULL,ControlTypeEnum);
     34                parameters->AddObject(new IntVecParam(ControlTypeEnum,control_type,num_control_type));
     35                xfree((void**)&control_type);
    3336
    3437                /*What solution type?*/
    … …  
    3942                        parameters->AddObject(new BoolParam(ControlSteadyEnum,false));
    4043                }
    41                 parameters->AddObject(new IntParam(NstepsEnum,iomodel->nsteps));
    42                 parameters->AddObject(new DoubleParam(TolxEnum,iomodel->tolx));
    43                 parameters->AddObject(new DoubleParam(EpsCmEnum,iomodel->eps_cm));
    44                 parameters->AddObject(new DoubleParam(MeanvelEnum,iomodel->meanvel));
    45 
    46                 parameters->AddObject(new BoolParam(CmGradientEnum,iomodel->cm_gradient));
    4744
    4845                /*Now, recover fit, optscal and maxiter as vectors: */
    49                 iomodel->FetchData(&iomodel->cm_responses,NULL,NULL,CmResponsesEnum);
    50                 iomodel->FetchData(&iomodel->cm_jump,NULL,NULL,CmJumpEnum);
     46                iomodel->FetchData(&iomodel->cm_responses,&nsteps,&num_control_type,CmResponsesEnum);
     47                iomodel->FetchData(&iomodel->cm_jump,&nsteps,&num_cm_responses,CmJumpEnum);
    5148                iomodel->FetchData(&iomodel->optscal,NULL,NULL,OptscalEnum);
    5249                iomodel->FetchData(&iomodel->maxiter,NULL,NULL,MaxiterEnum);
    5350
    54                 parameters->AddObject(new DoubleMatParam(OptscalEnum,iomodel->optscal,iomodel->nsteps,iomodel->num_control_type));
    55                 parameters->AddObject(new DoubleMatParam(CmResponsesEnum,iomodel->cm_responses,iomodel->nsteps,iomodel->num_cm_responses));
    56                 parameters->AddObject(new DoubleVecParam(CmJumpEnum,iomodel->cm_jump,iomodel->nsteps));
    57                 parameters->AddObject(new DoubleVecParam(MaxiterEnum,iomodel->maxiter,iomodel->nsteps));
     51                parameters->AddObject(new DoubleMatParam(OptscalEnum,iomodel->optscal,nsteps,num_control_type));
     52                parameters->AddObject(new DoubleMatParam(CmResponsesEnum,iomodel->cm_responses,nsteps,num_cm_responses));
     53                parameters->AddObject(new DoubleVecParam(CmJumpEnum,iomodel->cm_jump,nsteps));
     54                parameters->AddObject(new DoubleVecParam(MaxiterEnum,iomodel->maxiter,nsteps));
    5855
    5956                xfree((void**)&iomodel->cm_responses);

  • issm/trunk/src/c/modules/ModelProcessorx/Control/UpdateElementsAndMaterialsControl.cpp

    r9343 r9356  
    2121        Element  *element = NULL;
    2222        Material *material = NULL;
     23        int    numberofelements;
     24        int    num_control_type;
     25        bool   control_analysis;
     26
     27       
     28        /*Fetch parameters: */
     29        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     30        iomodel->parameters->FindParam(&control_analysis,ControlAnalysisEnum);
     31        iomodel->parameters->FindParam(&num_control_type,NumControlTypeEnum);
    2332
    2433        /*Now, return if no control*/
    25         if (!iomodel->control_analysis) return;
     34        if (!control_analysis) return;
    2635
    2736        /*Fetch data needed: */
    2837        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
    29 
    3038        iomodel->FetchDataToInput(elements,VxObsEnum);
    3139        iomodel->FetchDataToInput(elements,VyObsEnum);
    … …  
    3644        iomodel->FetchData(&iomodel->cm_min,NULL,NULL,CmMinEnum);
    3745        iomodel->FetchData(&iomodel->cm_max,NULL,NULL,CmMaxEnum);
    38         for(i=0;i<iomodel->num_control_type;i++){
     46        for(i=0;i<num_control_type;i++){
    3947                switch((int)iomodel->control_type[i]){
    4048                        case DhdtEnum:
    … …  
    6068        /*Update elements and materials: */
    6169        counter=0;
    62         for (i=0;i<iomodel->numberofelements;i++){
     70        for (i=0;i<numberofelements;i++){
    6371                if(iomodel->my_elements[i]){
    6472                        element=(Element*)elements->GetObjectByOffset(counter);

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

    r9340 r9356  
    1717        /*Intermediary*/
    1818        int i,j,k,n;
     19        int    dim;
     20        int    numberofelements;
     21        int    numberofvertices;
    1922
    2023        /*DataSets: */
    … …  
    2225        Vertices*     vertices = NULL;
    2326        Materials*    materials = NULL;
     27
     28        /*Fetch parameters: */
     29        iomodel->parameters->FindParam(&dim,DimEnum);
     30        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     31        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    2432
    2533        /*Did we already create the elements? : */
    … …  
    4654       
    4755        /*Create elements and materials: */
    48         for (i=0;i<iomodel->numberofelements;i++){
     56        for (i=0;i<numberofelements;i++){
    4957
    5058                if(iomodel->my_elements[i]){
    5159
    5260                        /*Create and add tria element to elements dataset: */
    53                         if(iomodel->dim==2)
     61                        if(dim==2)
    5462                         elements->AddObject(new Tria(i+1,i,i,iomodel,nummodels));
    5563                        else
    … …  
    7381
    7482        /*Add new constrant material property tgo materials, at the end: */
    75         materials->AddObject(new Matpar(iomodel->numberofelements+1,iomodel));//put it at the end of the materials
     83        materials->AddObject(new Matpar(numberofelements+1,iomodel));//put it at the end of the materials
    7684       
    7785        /*Create vertices: */
    … …  
    8290        iomodel->FetchData(&iomodel->thickness,NULL,NULL,ThicknessEnum);
    8391       
    84         for (i=0;i<iomodel->numberofvertices;i++){
     92        for (i=0;i<numberofvertices;i++){
    8593
    8694                /*vertices and nodes (same number, as we are running continuous galerkin formulation: */

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

    r9320 r9356  
    2020        int vertexid;
    2121        int elementswidth;
     22        int    dim;
     23        int    numberofelements;
     24        int    numberofvertices;
    2225
    2326        /*output*/
    2427        int* connectivity=NULL;
    2528
     29        /*Fetch parameters: */
     30        iomodel->parameters->FindParam(&dim,DimEnum);
     31        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     32        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     33
    2634        /*Some checks if debugging*/
    27         _assert_(iomodel->numberofvertices);
    28         _assert_(iomodel->numberofelements);
     35        _assert_(numberofvertices);
     36        _assert_(numberofelements);
    2937        _assert_(iomodel->elements);
    3038
    3139        /*Allocate ouput*/
    32         connectivity=(int*)xcalloc(iomodel->numberofvertices,sizeof(int));
     40        connectivity=(int*)xcalloc(numberofvertices,sizeof(int));
    3341
    3442        /*Get element width (3 or 6)*/
    35         if (iomodel->dim==2){
     43        if (dim==2){
    3644                elementswidth=3;
    3745        }
    … …  
    4149
    4250        /*Create connectivity table*/
    43         for (i=0;i<iomodel->numberofelements;i++){
     51        for (i=0;i<numberofelements;i++){
    4452                for (j=0;j<elementswidth;j++){
    4553                        vertexid=(int)iomodel->elements[elementswidth*i+j];
    46                         _assert_(vertexid>0 && vertexid-1<iomodel->numberofvertices);
     54                        _assert_(vertexid>0 && vertexid-1<numberofvertices);
    4755                        connectivity[vertexid-1]+=1;
    4856                }

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

    r9340 r9356  
    2424               
    2525        /*Initialize dataset: */
    26         parameters = new Parameters(ParametersEnum);
     26        parameters = (Parameters*)iomodel->parameters->Copy();
    2727
     28        /*some parameters that did not come with the iomodel: */
    2829        parameters->AddObject(new IntParam(SolutionTypeEnum,solution_type));
    2930        parameters->AddObject(new IntParam(AnalysisTypeEnum,analysis_type));
    3031        parameters->AddObject(new IntParam(AnalysisCounterEnum,analysis_counter));
    3132
    32         parameters->AddObject(new IntParam(DimEnum,iomodel->dim));
    33         parameters->AddObject(new BoolParam(IshutterEnum,iomodel->ishutter));
    34         parameters->AddObject(new BoolParam(IsmacayealpattynEnum,iomodel->ismacayealpattyn));
    35         parameters->AddObject(new BoolParam(IsstokesEnum,iomodel->isstokes));
    36         parameters->AddObject(new IntParam(OutputFrequencyEnum,iomodel->output_frequency));
    37         parameters->AddObject(new DoubleParam(EpsResEnum,iomodel->eps_res));
    38         parameters->AddObject(new DoubleParam(EpsRelEnum,iomodel->eps_rel));
    39         parameters->AddObject(new DoubleParam(EpsAbsEnum,iomodel->eps_abs));
    40         parameters->AddObject(new IntParam(MaxNonlinearIterationsEnum,(IssmInt)iomodel->max_nonlinear_iterations));
    41         parameters->AddObject(new IntParam(MaxSteadystateIterationsEnum,(IssmInt)iomodel->max_steadystate_iterations));
    42         parameters->AddObject(new DoubleParam(EpsvelEnum,iomodel->epsvel));
    43         parameters->AddObject(new DoubleParam(YtsEnum,iomodel->yts));
    44         parameters->AddObject(new DoubleParam(DtEnum,iomodel->dt*iomodel->yts));  //Ndt is in yeats in iomodel
    45         parameters->AddObject(new DoubleParam(NdtEnum,iomodel->ndt*iomodel->yts));//Ndt is in yeats in iomodel
    46         parameters->AddObject(new BoolParam(TimeAdaptEnum,iomodel->time_adapt));
    47         parameters->AddObject(new DoubleParam(CflCoefficientEnum,iomodel->cfl_coefficient));
    48         parameters->AddObject(new IntParam(HydrostaticAdjustmentEnum,iomodel->hydrostatic_adjustment));
    49         parameters->AddObject(new DoubleParam(PenaltyOffsetEnum,iomodel->penalty_offset));
    50         parameters->AddObject(new DoubleParam(SparsityEnum,iomodel->sparsity));
    51         parameters->AddObject(new BoolParam(LowmemEnum,iomodel->lowmem));
    52         parameters->AddObject(new IntParam(ConnectivityEnum,iomodel->connectivity));
    53         parameters->AddObject(new DoubleParam(BetaEnum,iomodel->beta));
    54         parameters->AddObject(new DoubleParam(MeltingpointEnum,iomodel->meltingpoint));
    55         parameters->AddObject(new DoubleParam(ReferencetemperatureEnum,iomodel->referencetemperature));
    56         parameters->AddObject(new DoubleParam(LatentheatEnum,iomodel->latentheat));
    57         parameters->AddObject(new DoubleParam(HeatcapacityEnum,iomodel->heatcapacity));
    58         parameters->AddObject(new IntParam(ArtDiffEnum,iomodel->artdiff));
    59         parameters->AddObject(new IntParam(BasalMeltingRateCorrectionApplyEnum,iomodel->basal_melting_rate_correction_apply));
    60         parameters->AddObject(new DoubleParam(GroundingLineMeltingRateEnum,iomodel->gl_melting_rate));
    61         parameters->AddObject(new DoubleParam(PenaltyMeltingEnum,iomodel->penalty_melting));
    62         parameters->AddObject(new IntParam(MinThermalConstraintsEnum,iomodel->min_thermal_constraints));
    63         parameters->AddObject(new IntParam(MinMechanicalConstraintsEnum,iomodel->min_mechanical_constraints));
    64         parameters->AddObject(new IntParam(StabilizeConstraintsEnum,iomodel->stabilize_constraints));
    65         parameters->AddObject(new DoubleParam(StokesreconditioningEnum,iomodel->stokesreconditioning));
    66         parameters->AddObject(new IntParam(ShelfDampeningEnum,iomodel->shelf_dampening));
    67         parameters->AddObject(new DoubleParam(ViscosityOvershootEnum,iomodel->viscosity_overshoot));
    68         parameters->AddObject(new BoolParam(WaitonlockEnum,iomodel->waitonlock));
    69         parameters->AddObject(new IntParam(NumberOfElementsEnum,iomodel->numberofelements));
    70         parameters->AddObject(new BoolParam(IoGatherEnum,iomodel->io_gather));
    71         parameters->AddObject(new IntParam(GroundingLineMigrationEnum,iomodel->gl_migration));
    72         parameters->AddObject(new IntParam(RheologyLawEnum,iomodel->rheology_law));
    73         parameters->AddObject(new BoolParam(IsdiagnosticEnum,iomodel->isdiagnostic));
    74         parameters->AddObject(new BoolParam(IsprognosticEnum,iomodel->isprognostic));
    75         parameters->AddObject(new BoolParam(IsthermalEnum,iomodel->isthermal));
    7633        parameters->AddObject(new DoubleParam(TimeEnum,0.0));  //start at time 0 by default for all solutions.
    77 
    78         /*Deal with more complex parameters*/
    79         iomodel->FetchData(&iomodel->riftinfo,&iomodel->numrifts,NULL,RiftinfoEnum);
    80         parameters->AddObject(new IntParam(NumRiftsEnum,iomodel->numrifts));
    81         xfree((void**)&iomodel->riftinfo);
    82 
    8334
    8435        /*Requested output?*/
    … …  
    9243        CreateParametersQmu(&parameters,iomodel,solution_type,analysis_type);
    9344
     45        /*Go through all parameters, and convert units to SI: */
     46        parameters->UnitConversion(ExtToIuEnum);
     47
    9448        /*Assign output pointer: */
    9549        *pparameters=parameters;

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

    r9320 r9356  
    2020        int vertexid;
    2121        int elementswidth;
     22        int    dim;
     23        int    numberofelements;
     24        int    numberofvertices;
    2225
    2326        /*output*/
    2427        int* connectivity=NULL;
    2528
     29        /*Fetch parameters: */
     30        iomodel->parameters->FindParam(&dim,DimEnum);
     31        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     32        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     33
    2634        /*Some checks if debugging*/
    27         _assert_(iomodel->numberofvertices);
    28         _assert_(iomodel->numberofelements);
     35        _assert_(numberofvertices);
     36        _assert_(numberofelements);
    2937        _assert_(iomodel->my_elements);
    3038        _assert_(iomodel->elements);
    3139
    3240        /*Allocate ouput*/
    33         connectivity=(int*)xcalloc(iomodel->numberofvertices,sizeof(int));
     41        connectivity=(int*)xcalloc(numberofvertices,sizeof(int));
    3442
    3543        /*Get element width (3 or 6)*/
    36         if (iomodel->dim==2){
     44        if (dim==2){
    3745                elementswidth=3;
    3846        }
    … …  
    4250
    4351        /*Create connectivity table*/
    44         for (i=0;i<iomodel->numberofelements;i++){
     52        for (i=0;i<numberofelements;i++){
    4553                /*!! in parallel we do not want the vertex to be connected to an element that is not in its partition!!*/
    4654                if(iomodel->my_elements[i]){
    4755                        for (j=0;j<elementswidth;j++){
    4856                                vertexid=(int)iomodel->elements[elementswidth*i+j];
    49                                 _assert_(vertexid>0 && vertexid-1<iomodel->numberofvertices);
     57                                _assert_(vertexid>0 && vertexid-1<numberofvertices);
    5058                                connectivity[vertexid-1]=i+1;
    5159                        }

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp

    r9340 r9356  
    1616        int i,j;
    1717        int count;
     18        double yts;
     19        double g;
     20        double rho_ice;
     21        double stokesreconditioning;
     22        bool   isstokes,ismacayealpattyn;
    1823       
    1924        /*Output*/
    … …  
    2126        SpcStatic*    spcstatic  = NULL;
    2227        int     node1,node2;
     28        int    dim;
     29        int    numberofvertices;
     30
     31        /*Fetch parameters: */
     32        iomodel->parameters->FindParam(&dim,DimEnum);
     33        iomodel->parameters->FindParam(&yts,YtsEnum);
     34        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     35        iomodel->parameters->FindParam(&g,GEnum);
     36        iomodel->parameters->FindParam(&rho_ice,RhoIceEnum);
     37        iomodel->parameters->FindParam(&stokesreconditioning,StokesreconditioningEnum);
     38        iomodel->parameters->FindParam(&isstokes,IsstokesEnum);
     39        iomodel->parameters->FindParam(&ismacayealpattyn,IsmacayealpattynEnum);
    2340
    2441        /*Recover pointer: */
    … …  
    2946       
    3047        /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */
    31         if (!iomodel->ismacayealpattyn & !iomodel->isstokes)goto cleanup_and_return;
     48        if (!ismacayealpattyn & !isstokes)goto cleanup_and_return;
    3249       
    3350        /*Constraints: fetch data: */
    … …  
    3754        iomodel->FetchData(&iomodel->nodeonhutter,NULL,NULL,NodeOnHutterEnum);
    3855        iomodel->FetchData(&iomodel->nodeonmacayeal,NULL,NULL,NodeOnMacayealEnum);
    39         if(iomodel->dim==3)iomodel->FetchData(&iomodel->nodeonpattyn,NULL,NULL,NodeOnPattynEnum);
    40         if(iomodel->dim==3)iomodel->FetchData(&iomodel->nodeonstokes,NULL,NULL,NodeOnStokesEnum);
     56        if(dim==3)iomodel->FetchData(&iomodel->nodeonpattyn,NULL,NULL,NodeOnPattynEnum);
     57        if(dim==3)iomodel->FetchData(&iomodel->nodeonstokes,NULL,NULL,NodeOnStokesEnum);
    4158        iomodel->FetchData(&iomodel->vertices_type,NULL,NULL,VerticesTypeEnum);
    4259        iomodel->FetchData(&iomodel->surface,NULL,NULL,SurfaceEnum);
    … …  
    4764       
    4865        /*Create spcs from x,y,z, as well as the spc values on those spcs: */
    49         for (i=0;i<iomodel->numberofvertices;i++){
     66        for (i=0;i<numberofvertices;i++){
    5067                if(iomodel->my_vertices[i]){
    5168
    … …  
    5976                                                count++;
    6077                                                if (!isnan(iomodel->spcvx[i])){
    61                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    62                                                         count++;
    63                                                 }
    64                                                 if (!isnan(iomodel->spcvy[i])){
    65                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     78                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     79                                                        count++;
     80                                                }
     81                                                if (!isnan(iomodel->spcvy[i])){
     82                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    6683                                                        count++;
    6784                                                }
    … …  
    7491                                                count++;
    7592                                                if (!isnan(iomodel->spcvx[i])){
    76                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    77                                                         count++;
    78                                                 }
    79                                                 if (!isnan(iomodel->spcvy[i])){
    80                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     93                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     94                                                        count++;
     95                                                }
     96                                                if (!isnan(iomodel->spcvy[i])){
     97                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    8198                                                        count++;
    8299                                                }
    … …  
    96113                                                count++;
    97114                                                if (!isnan(iomodel->spcvx[i])){
    98                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    99                                                         count++;
    100                                                 }
    101                                                 if (!isnan(iomodel->spcvy[i])){
    102                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     115                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     116                                                        count++;
     117                                                }
     118                                                if (!isnan(iomodel->spcvy[i])){
     119                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    103120                                                        count++;
    104121                                                }
    … …  
    111128                                                count++;
    112129                                                if (!isnan(iomodel->spcvx[i])){
    113                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    114                                                         count++;
    115                                                 }
    116                                                 if (!isnan(iomodel->spcvy[i])){
    117                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     130                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     131                                                        count++;
     132                                                }
     133                                                if (!isnan(iomodel->spcvy[i])){
     134                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    118135                                                        count++;
    119136                                                }
    120137                                                if (!isnan(iomodel->spcvz[i])){
    121                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,iomodel->spcvz[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     138                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,iomodel->spcvz[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    122139                                                        count++;
    123140                                                }
    … …  
    136153                                                count++;
    137154                                                if (!isnan(iomodel->spcvx[i])){
    138                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    139                                                         count++;
    140                                                 }
    141                                                 if (!isnan(iomodel->spcvy[i])){
    142                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     155                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     156                                                        count++;
     157                                                }
     158                                                if (!isnan(iomodel->spcvy[i])){
     159                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    143160                                                        count++;
    144161                                                }
    … …  
    151168                                                count++;
    152169                                                if (!isnan(iomodel->spcvx[i])){
    153                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    154                                                         count++;
    155                                                 }
    156                                                 if (!isnan(iomodel->spcvy[i])){
    157                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     170                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     171                                                        count++;
     172                                                }
     173                                                if (!isnan(iomodel->spcvy[i])){
     174                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    158175                                                        count++;
    159176                                                }
    160177                                                if (!isnan(iomodel->spcvz[i])){
    161                                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,iomodel->spcvz[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     178                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,5,iomodel->spcvz[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    162179                                                        count++;
    163180                                                }
    … …  
    168185                        else{
    169186                                if (!isnan(iomodel->spcvx[i])){
    170                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     187                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    171188                                        count++;
    172189                                }
    … …  
    177194                               
    178195                                if (!isnan(iomodel->spcvy[i])){
    179                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vy.
     196                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vy.
    180197                                        count++;
    181198                                }
    … …  
    186203
    187204                                if (!isnan(iomodel->spcvz[i]) && ((int)iomodel->vertices_type[i]==StokesApproximationEnum ||  ((int)iomodel->vertices_type[i]==NoneApproximationEnum))){
    188                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvz[i]/iomodel->yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy
     205                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,3,iomodel->spcvz[i]/yts,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy
    189206                                        count++;
    190207                                }
    191208                                if ((int)iomodel->vertices_type[i]==NoneApproximationEnum){
    192                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,iomodel->g*iomodel->rho_ice*(iomodel->surface[i]-iomodel->z[i])/iomodel->stokesreconditioning,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy
     209                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,4,g*rho_ice*(iomodel->surface[i]-iomodel->z[i])/stokesreconditioning,DiagnosticHorizAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy
    193210                                        count++;
    194211                                }

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateLoadsDiagnosticHoriz.cpp

    r9340 r9356  
    2525        int count=0;
    2626        int penpair_ids[2];
     27        int dim;
     28        int numberofvertices;
     29        bool ismacayealpattyn,isstokes;
     30        int  numpenalties,numberofpressureloads,numrifts;
     31
     32        /*Fetch parameters: */
     33        iomodel->parameters->FindParam(&dim,DimEnum);
     34        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     35        iomodel->parameters->FindParam(&isstokes,IsstokesEnum);
     36        iomodel->parameters->FindParam(&ismacayealpattyn,IsmacayealpattynEnum);
    2737
    2838        /*Recover pointer: */
    … …  
    3343
    3444        /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */
    35         if (!iomodel->ismacayealpattyn & !iomodel->isstokes)goto cleanup_and_return;
     45        if (!ismacayealpattyn & !isstokes)goto cleanup_and_return;
    3646       
    3747        /*Create pressure loads as boundary conditions. Pay attention to the partitioning if we are running in parallel (the nodes
    3848         * referenced by a certain load must belong to the cluster node): */
    39         iomodel->FetchData(&iomodel->pressureload,&iomodel->numberofpressureloads,NULL,PressureloadEnum);
     49        iomodel->FetchData(&iomodel->pressureload,&numberofpressureloads,NULL,PressureloadEnum);
    4050        iomodel->FetchData(&iomodel->elements_type,NULL,NULL,ElementsTypeEnum);
    4151        iomodel->FetchData(&iomodel->thickness,NULL,NULL,ThicknessEnum);
    … …  
    4656
    4757        /*First load data:*/
    48         for (i=0;i<iomodel->numberofpressureloads;i++){
     58        for (i=0;i<numberofpressureloads;i++){
    4959               
    5060                /*Retrieve element to which this icefront belongs: */
    51                 if (iomodel->dim==2) segment_width=4;
     61                if (dim==2) segment_width=4;
    5262                else segment_width=6;
    5363                element=(int)(*(iomodel->pressureload+segment_width*i+segment_width-2)-1); //element is in the penultimate column (node1 node2 ... elem fill)
    … …  
    6070
    6171                /*Create and  add load: */
    62                 if ((int)*(iomodel->elements_type+element)==(MacAyealApproximationEnum) && iomodel->dim==2){
     72                if ((int)*(iomodel->elements_type+element)==(MacAyealApproximationEnum) && dim==2){
    6373                        loads->AddObject(new Icefront(iomodel->loadcounter+count+1,i,iomodel,MacAyeal2dIceFrontEnum,DiagnosticHorizAnalysisEnum));
    6474                        count++;
    6575                }
    66                 else if ((int)*(iomodel->elements_type+element)==(MacAyealApproximationEnum) && iomodel->dim==3){
     76                else if ((int)*(iomodel->elements_type+element)==(MacAyealApproximationEnum) && dim==3){
    6777                        loads->AddObject(new Icefront(iomodel->loadcounter+count+1,i,iomodel,MacAyeal3dIceFrontEnum,DiagnosticHorizAnalysisEnum));
    6878                        count++;
    … …  
    113123        CreateSingleNodeToElementConnectivity(iomodel);
    114124       
    115         for (i=0;i<iomodel->numberofvertices;i++){
     125        for (i=0;i<numberofvertices;i++){
    116126               
    117127                if(iomodel->my_vertices[i]==1 && iomodel->singlenodetoelementconnectivity[i]!=0){
    … …  
    137147
    138148        /*Create Penpair for penalties: */
    139         iomodel->FetchData(&iomodel->penalties,&iomodel->numpenalties,NULL,PenaltiesEnum);
     149        iomodel->FetchData(&iomodel->penalties,&numpenalties,NULL,PenaltiesEnum);
    140150       
    141         for(i=0;i<iomodel->numpenalties;i++){
     151        for(i=0;i<numpenalties;i++){
    142152
    143153                if(iomodel->my_vertices[(int)iomodel->penalties[2*i+0]-1]){
    … …  
    160170
    161171        /*Create Riffront loads for rifts: */
    162         iomodel->FetchData(&iomodel->riftinfo,&iomodel->numrifts,NULL,RiftinfoEnum);
     172        iomodel->FetchData(&iomodel->riftinfo,&numrifts,NULL,RiftinfoEnum);
    163173        iomodel->FetchData(&iomodel->thickness,NULL,NULL,ThicknessEnum);
    164174        iomodel->FetchData(&iomodel->bed,NULL,NULL,BedEnum);
    … …  
    167177
    168178       
    169         for(i=0;i<iomodel->numrifts;i++){
     179        for(i=0;i<numrifts;i++){
    170180
    171181                if(iomodel->my_elements[(int)*(iomodel->riftinfo+RIFTINFOSIZE*i+2)-1]){

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateNodesDiagnosticHoriz.cpp

    r9340 r9356  
    1717        /*Intermediary*/
    1818        int i;
    19         bool continuous_galerkin=true;
     19        bool   continuous_galerkin=true;
     20        int    numberofvertices;
     21        bool   isstokes,ismacayealpattyn;
    2022
    2123        /*DataSets: */
    2224        Nodes*    nodes = NULL;
     25
     26        /*Fetch parameters: */
     27        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     28        iomodel->parameters->FindParam(&isstokes,IsstokesEnum);
     29        iomodel->parameters->FindParam(&ismacayealpattyn,IsmacayealpattynEnum);
    2330
    2431        /*Recover pointer: */
    … …  
    2936       
    3037        /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */
    31         if (!iomodel->ismacayealpattyn & !iomodel->isstokes)goto cleanup_and_return;
     38        if (!ismacayealpattyn & !isstokes)goto cleanup_and_return;
    3239
    3340        /*Continuous Galerkin partition of nodes: */
    … …  
    4653        iomodel->FetchData(&iomodel->diagnostic_ref,NULL,NULL,DiagnosticRefEnum);
    4754       
    48         for (i=0;i<iomodel->numberofvertices;i++){
     55        for (i=0;i<numberofvertices;i++){
    4956
    5057                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHoriz/UpdateElementsDiagnosticHoriz.cpp

    r9343 r9356  
    1616void    UpdateElementsDiagnosticHoriz(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
    18         /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */
    19         if (!iomodel->ismacayealpattyn & !iomodel->isstokes) return;
    20 
     18        int    dim;
     19        int    numberofelements;
     20        bool   ismacayealpattyn;
     21        bool   isstokes;
     22        bool   control_analysis;
     23        bool   qmu_analysis;
     24               
    2125        /*Fetch data needed: */
    2226        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
    2327        iomodel->FetchData(&iomodel->elements_type,NULL,NULL,ElementsTypeEnum);
    2428
     29        iomodel->parameters->FindParam(&isstokes,IsstokesEnum);
     30        iomodel->parameters->FindParam(&ismacayealpattyn,IsmacayealpattynEnum);
     31        iomodel->parameters->FindParam(&dim,DimEnum);
     32        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     33        iomodel->parameters->FindParam(&control_analysis,ControlAnalysisEnum);
     34        iomodel->parameters->FindParam(&qmu_analysis,QmuAnalysisEnum);
     35
     36
     37        /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */
     38        if (!ismacayealpattyn & !isstokes) return;
     39
    2540        /*Update elements: */
    2641        int counter=0;
    27         for(int i=0;i<iomodel->numberofelements;i++){
     42        for(int i=0;i<numberofelements;i++){
    2843                if(iomodel->my_elements[i]){
    2944                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    5065        iomodel->FetchDataToInput(elements,BathymetryEnum);
    5166
    52         if (iomodel->dim==3){
     67        if (dim==3){
    5368                iomodel->FetchDataToInput(elements,ElementOnBedEnum);
    5469                iomodel->FetchDataToInput(elements,ElementOnSurfaceEnum);
    … …  
    6176        }
    6277
    63         if(iomodel->control_analysis){
     78        if(control_analysis){
    6479                iomodel->FetchDataToInput(elements,VxObsEnum);
    6580                iomodel->FetchDataToInput(elements,VyObsEnum);
    … …  
    7085        elements->InputDuplicate(VxEnum,VxPicardEnum);
    7186        elements->InputDuplicate(VxEnum,VxObsEnum);
    72         if(iomodel->qmu_analysis)elements->InputDuplicate(VxEnum,QmuVxEnum);
     87        if(qmu_analysis)elements->InputDuplicate(VxEnum,QmuVxEnum);
    7388       
    7489        elements->InputDuplicate(VyEnum,VyPicardEnum);
    7590        elements->InputDuplicate(VyEnum,VyObsEnum);
    76         if(iomodel->qmu_analysis)elements->InputDuplicate(VyEnum,QmuVyEnum);
     91        if(qmu_analysis)elements->InputDuplicate(VyEnum,QmuVyEnum);
    7792       
    78         if(iomodel->dim==3){
     93        if(dim==3){
    7994                elements->InputDuplicate(VzEnum,VzPicardEnum);
    8095                elements->InputDuplicate(VzEnum,VzObsEnum);
    81                 if(iomodel->qmu_analysis)elements->InputDuplicate(VzEnum,QmuVzEnum);
     96                if(qmu_analysis)elements->InputDuplicate(VzEnum,QmuVzEnum);
    8297        }
    8398       

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateConstraintsDiagnosticHutter.cpp

    r9340 r9356  
    1616        int i;
    1717        int count;
     18        double yts;
     19        int    numberofvertices;
     20        bool   ishutter;
    1821
    1922        /*Output*/
    … …  
    2427        constraints=*pconstraints;
    2528
     29        /*Fetch parameters: */
     30        iomodel->parameters->FindParam(&yts,YtsEnum);
     31        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     32        iomodel->parameters->FindParam(&ishutter,IshutterEnum);
     33
    2634        /*Create constraints if they do not exist yet*/
    2735        if(!constraints) constraints = new Constraints(ConstraintsEnum);
    2836
    2937        /*Now, is the flag ishutter on? otherwise, do nothing: */
    30         if (!iomodel->ishutter) goto cleanup_and_return;
     38        if (!ishutter) goto cleanup_and_return;
    3139
    3240        /*Fetch data: */
    … …  
    3947
    4048        /*vx and vy are spc'd if we are not on nodeonhutter: */
    41         for (i=0;i<iomodel->numberofvertices;i++){
     49        for (i=0;i<numberofvertices;i++){
    4250                /*keep only this partition's nodes:*/
    4351                if((iomodel->my_vertices[i])){
    … …  
    5260                        else{
    5361                                if (!isnan(iomodel->spcvx[i])){
    54                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/iomodel->yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     62                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->spcvx[i]/yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    5563                                        count++;
    5664                                }
    5765
    5866                                if (!isnan(iomodel->spcvy[i])){
    59                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/iomodel->yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy
     67                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,2,iomodel->spcvy[i]/yts,DiagnosticHutterAnalysisEnum)); //add count'th spc, on node i+1, setting dof 2 to vy
    6068                                        count++;
    6169                                }

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/CreateNodesDiagnosticHutter.cpp

    r9340 r9356  
    1818        int i;
    1919        bool continuous_galerkin=true;
     20        int    numberofvertices;
     21        bool   ishutter;
    2022
    2123        /*DataSets: */
    2224        Nodes*    nodes = NULL;
     25
     26        /*Fetch parameters: */
     27        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     28        iomodel->parameters->FindParam(&ishutter,IshutterEnum);
    2329
    2430        /*Recover pointer: */
    … …  
    2935
    3036        /*Now, is the flag ishutter on? otherwise, do nothing: */
    31         if (!iomodel->ishutter)goto cleanup_and_return;
     37        if (!ishutter)goto cleanup_and_return;
    3238
    3339        /*Continuous Galerkin partition of nodes: */
    … …  
    4652        CreateNumberNodeToElementConnectivity(iomodel);
    4753
    48         for (i=0;i<iomodel->numberofvertices;i++){
     54        for (i=0;i<numberofvertices;i++){
    4955
    5056                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticHutter/UpdateElementsDiagnosticHutter.cpp

    r9343 r9356  
    1616void    UpdateElementsDiagnosticHutter(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
    18         /*Now, is the flag hutter on? otherwise, do nothing: */
    19         if (!iomodel->ishutter)return;
     18        int    numberofelements;
     19        bool   ishutter;
    2020
     21       
    2122        /*Fetch data needed: */
    2223        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
    2324        iomodel->FetchData(&iomodel->elements_type,NULL,NULL,ElementsTypeEnum);
     25        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     26        iomodel->parameters->FindParam(&ishutter,IshutterEnum);
     27
     28        /*Now, is the flag hutter on? otherwise, do nothing: */
     29        if (!ishutter)return;
    2430
    2531        /*Update elements: */
    2632        int counter=0;
    27         for(int i=0;i<iomodel->numberofelements;i++){
     33        for(int i=0;i<numberofelements;i++){
    2834                if(iomodel->my_elements[i]){
    2935                        Element* element=(Element*)elements->GetObjectByOffset(counter);

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateConstraintsDiagnosticVert.cpp

    r9340 r9356  
    1515        /*Intermediary*/
    1616        int i;
     17        int dim;
    1718        int count;
     19        double yts;
     20        int    numberofvertices;
    1821
    1922        /*Output*/
    2023        Constraints* constraints = NULL;
     24
     25        /*Fetch parameters: */
     26        iomodel->parameters->FindParam(&dim,DimEnum);
     27        iomodel->parameters->FindParam(&yts,YtsEnum);
    2128
    2229        /*Recover pointer: */
    … …  
    2734
    2835        /*return if 2d mesh*/
    29         if (iomodel->dim==2) goto cleanup_and_return;
     36        if (dim==2) goto cleanup_and_return;
    3037
    3138        /*Fetch data: */
    3239        iomodel->FetchData(&iomodel->spcvz,NULL,NULL,SpcvzEnum);
    3340        iomodel->FetchData(&iomodel->nodeonstokes,NULL,NULL,NodeOnStokesEnum);
     41        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    3442
    3543        /*Initialize counter*/
    … …  
    3745
    3846        /*Create spcs from x,y,z, as well as the spc values on those spcs: */
    39         for (i=0;i<iomodel->numberofvertices;i++){
     47        for (i=0;i<numberofvertices;i++){
    4048
    4149                /*keep only this partition's nodes:*/
    … …  
    4856                        else if (!isnan(iomodel->spcvz[i])){
    4957                                constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,
    50                                                                 iomodel->spcvz[i]/iomodel->yts,DiagnosticVertAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
     58                                                                iomodel->spcvz[i]/yts,DiagnosticVertAnalysisEnum)); //add count'th spc, on node i+1, setting dof 1 to vx.
    5159                                count++;
    5260

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/CreateNodesDiagnosticVert.cpp

    r9340 r9356  
    1818        int i;
    1919        bool continuous_galerkin=true;
     20        int    dim;
     21        int    numberofvertices;
    2022
    2123        /*DataSets: */
    2224        Nodes*    nodes = NULL;
     25
     26        /*Fetch parameters: */
     27        iomodel->parameters->FindParam(&dim,DimEnum);
     28        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    2329
    2430        /*Recover pointer: */
    … …  
    2935
    3036        /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */
    31         if (iomodel->dim==2)goto cleanup_and_return;
     37        if (dim==2)goto cleanup_and_return;
    3238
    3339        /*Continuous Galerkin partition of nodes: */
    … …  
    4248        iomodel->FetchData(&iomodel->nodeonwater,NULL,NULL,NodeOnWaterEnum);
    4349
    44         for (i=0;i<iomodel->numberofvertices;i++){
     50        for (i=0;i<numberofvertices;i++){
    4551
    4652                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/DiagnosticVert/UpdateElementsDiagnosticVert.cpp

    r9343 r9356  
    1616void    UpdateElementsDiagnosticVert(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    dim;
     19        int    numberofelements;
     20
     21        /*Fetch parameters: */
     22        iomodel->parameters->FindParam(&dim,DimEnum);
     23        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     24
    1825        /*Now, is the model 3d? otherwise, do nothing: */
    19         if (iomodel->dim==2)return;
     26        if (dim==2)return;
    2027
    2128        /*Fetch data needed: */
    … …  
    2431        /*Update elements: */
    2532        int counter=0;
    26         for(int i=0;i<iomodel->numberofelements;i++){
     33        for(int i=0;i<numberofelements;i++){
    2734                if(iomodel->my_elements[i]){
    2835                        Element* element=(Element*)elements->GetObjectByOffset(counter);

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

    r9340 r9356  
    2424        extern int my_rank;
    2525        extern int num_procs;
     26        int    numberofelements;
     27        int    numberofvertices;
     28        int    numberofelements2d;
     29        int    numberofvertices2d;
     30        int    numlayers;
     31        int    numrifts;
     32        int    numpenalties;
    2633
    2734        /*output: */
    … …  
    3441        int  elements_width; //number of columns in elements (2d->3, 3d->6)
    3542        int  el1,el2;
     43        int    dim;
     44
     45        /*Fetch parameters: */
     46        iomodel->parameters->FindParam(&dim,DimEnum);
     47        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     48        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     49        iomodel->parameters->FindParam(&numberofelements2d,NumberOfElements2DEnum);
     50        iomodel->parameters->FindParam(&numberofvertices2d,NumberOfNodes2DEnum);
     51        iomodel->parameters->FindParam(&numlayers,NumlayersEnum);
     52        iomodel->parameters->FindParam(&numrifts,NumRiftsEnum);
    3653
    3754        /*First, check that partitioning has not yet been carryed out. Just check whether my_elements pointers is not already assigned a value: */
    … …  
    3956
    4057        /*Number of vertices per elements, needed to correctly retrieve data: */
    41         if(iomodel->dim==2) elements_width=3; //tria elements
     58        if(dim==2) elements_width=3; //tria elements
    4259        else elements_width=6; //penta elements
    4360
    4461        #ifdef _PARALLEL_
    4562        /*Determine parallel partitioning of elements: we use Metis for now. First load the data, then partition*/
    46         if(iomodel->dim==2){
     63        if(dim==2){
    4764                /*load elements: */
    4865                iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
    … …  
    5370        }
    5471
    55         MeshPartitionx(&epart, &npart,iomodel->numberofelements,iomodel->numberofvertices,iomodel->elements, iomodel->numberofelements2d,iomodel->numberofvertices2d,iomodel->elements2d,iomodel->numlayers,elements_width, iomodel->dim,num_procs);
     72        MeshPartitionx(&epart, &npart,numberofelements,numberofvertices,iomodel->elements, numberofelements2d,numberofvertices2d,iomodel->elements2d,numlayers,elements_width, dim,num_procs);
    5673
    5774        /*Free elements and elements2d: */
    … …  
    6178        #else
    6279        /*In serial mode, epart is full of 0: all elements belong to cpu 0: */
    63         epart=(int*)xcalloc(iomodel->numberofelements,sizeof(int));
     80        epart=(int*)xcalloc(numberofelements,sizeof(int));
    6481        #endif
    6582
    6683        /*Deal with rifts, they have to be included into one partition only, not several: */
    67         iomodel->FetchData(&iomodel->riftinfo,&iomodel->numrifts,NULL,RiftinfoEnum);
     84        iomodel->FetchData(&iomodel->riftinfo,&numrifts,NULL,RiftinfoEnum);
    6885
    69         for(i=0;i<iomodel->numrifts;i++){
     86        for(i=0;i<numrifts;i++){
    7087                el1=(int)*(iomodel->riftinfo+RIFTINFOSIZE*i+2)-1; //matlab indexing to c indexing
    7188                el2=(int)*(iomodel->riftinfo+RIFTINFOSIZE*i+3)-1; //matlab indexing to c indexing
    … …  
    7794
    7895        /*Used later on: */
    79         my_vertices=(int*)xcalloc(iomodel->numberofvertices,sizeof(int));
    80         my_elements=(bool*)xcalloc(iomodel->numberofelements,sizeof(bool));
     96        my_vertices=(int*)xcalloc(numberofvertices,sizeof(int));
     97        my_elements=(bool*)xcalloc(numberofelements,sizeof(bool));
    8198
    8299        /*Start figuring out, out of the partition, which elements belong to this cpu: */
    83100        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
    84         for (i=0;i<iomodel->numberofelements;i++){
     101        for (i=0;i<numberofelements;i++){
    85102
    86103                /*!All elements have been partitioned above, only deal with elements for this cpu: */
    … …  
    110127         * penpair has 2 nodes that are poointing toward 2 vertices.
    111128         * The 2 vertices must be in the same cpu as the penpair*/
    112         iomodel->FetchData(&iomodel->penalties,&iomodel->numpenalties,NULL,PenaltiesEnum);
    113         for(i=0;i<iomodel->numpenalties;i++){
     129        iomodel->FetchData(&iomodel->penalties,&numpenalties,NULL,PenaltiesEnum);
     130        for(i=0;i<numpenalties;i++){
    114131                if(my_vertices[(int)iomodel->penalties[2*i+0]-1] && !my_vertices[(int)iomodel->penalties[2*i+1]-1]){
    115132                        my_vertices[(int)iomodel->penalties[2*i+1]-1]=2; //to know that these elements are not on the partition

  • issm/trunk/src/c/modules/ModelProcessorx/Enthalpy/CreateConstraintsEnthalpy.cpp

    r9340 r9356  
    1616        int i;
    1717        int count;
     18        int    dim;
     19        int    numberofvertices;
     20        double heatcapacity;
     21        double referencetemperature;
    1822       
    1923        /*Output*/
    2024        Constraints* constraints = NULL;
    2125        SpcStatic*    spcstatic  = NULL;
     26
     27        /*Fetch parameters: */
     28        iomodel->parameters->FindParam(&dim,DimEnum);
     29        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     30        iomodel->parameters->FindParam(&heatcapacity,HeatcapacityEnum);
     31        iomodel->parameters->FindParam(&referencetemperature,ReferencetemperatureEnum);
    2232
    2333        /*Recover pointer: */
    … …  
    2838
    2939        /*return if 2d mesh*/
    30         if (iomodel->dim==2) goto cleanup_and_return;
     40        if (dim==2) goto cleanup_and_return;
    3141
    3242        /*Fetch data: */
    … …  
    3747
    3848        /*Create constraints from x,y,z: */
    39         for (i=0;i<iomodel->numberofvertices;i++){
     49        for (i=0;i<numberofvertices;i++){
    4050                /*keep only this partition's nodes:*/
    4151                if((iomodel->my_vertices[i])){
    … …  
    4353                        if ((int)iomodel->spctemperature[2*i]){
    4454
    45                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,iomodel->heatcapacity*(iomodel->spctemperature[2*i+1]-iomodel->referencetemperature),EnthalpyAnalysisEnum));
     55                                constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,1,heatcapacity*(iomodel->spctemperature[2*i+1]-referencetemperature),EnthalpyAnalysisEnum));
    4656                                count++;
    4757

  • issm/trunk/src/c/modules/ModelProcessorx/Enthalpy/CreateNodesEnthalpy.cpp

    r9340 r9356  
    1818        int i;
    1919        bool continuous_galerkin=true;
     20        int    numberofvertices;
    2021
    2122        /*DataSets: */
    2223        Nodes*    nodes = NULL;
     24
     25        /*Fetch parameters: */
     26        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    2327       
    2428        /*Recover pointer: */
    … …  
    3943        iomodel->FetchData(&iomodel->nodeonwater,NULL,NULL,NodeOnWaterEnum);
    4044
    41         for (i=0;i<iomodel->numberofvertices;i++){
     45        for (i=0;i<numberofvertices;i++){
    4246
    4347                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/Enthalpy/UpdateElementsEnthalpy.cpp

    r9343 r9356  
    1616void    UpdateElementsEnthalpy(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    dim;
     19        int    numberofelements;
     20
     21        /*Fetch parameters: */
     22        iomodel->parameters->FindParam(&dim,DimEnum);
     23        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     24
    1825        /*Now, is the model 3d? otherwise, do nothing: */
    19         if (iomodel->dim==2)return;
     26        if (dim==2)return;
    2027
    2128        /*Fetch data needed: */
    … …  
    2431        /*Update elements: */
    2532        int counter=0;
    26         for(int i=0;i<iomodel->numberofelements;i++){
     33        for(int i=0;i<numberofelements;i++){
    2734                if(iomodel->my_elements[i]){
    2835                        Element* element=(Element*)elements->GetObjectByOffset(counter);

  • issm/trunk/src/c/modules/ModelProcessorx/Hydrology/CreateNodesHydrology.cpp

    r9340 r9356  
    1818        int i;
    1919        bool continuous_galerkin=true;
     20        int    numberofvertices;
    2021
    2122        /*DataSets: */
    2223        Nodes*    nodes = NULL;
     24
     25        /*Fetch parameters: */
     26        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    2327       
    2428        /*Recover pointer: */
    … …  
    3943        iomodel->FetchData(&iomodel->nodeonwater,NULL,NULL,NodeOnWaterEnum);
    4044
    41         for (i=0;i<iomodel->numberofvertices;i++){
     45        for (i=0;i<numberofvertices;i++){
    4246
    4347                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/Hydrology/UpdateElementsHydrology.cpp

    r9343 r9356  
    1616void    UpdateElementsHydrology(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    numberofelements;
     19       
    1820        /*Fetch data needed: */
    1921        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
     22        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
    2023
    2124        /*Update elements: */
    2225        int counter=0;
    23         for(int i=0;i<iomodel->numberofelements;i++){
     26        for(int i=0;i<numberofelements;i++){
    2427                if(iomodel->my_elements[i]){
    2528                        Element* element=(Element*)elements->GetObjectByOffset(counter);

  • issm/trunk/src/c/modules/ModelProcessorx/Melting/CreateLoadsMelting.cpp

    r9340 r9356  
    1515        /*Intermediary*/
    1616        int i;
     17        int    dim;
     18        int    numberofvertices;
    1719
    1820        /*DataSet*/
    1921        Loads*    loads    = NULL;
    2022
     23        /*Fetch parameters: */
     24        iomodel->parameters->FindParam(&dim,DimEnum);
     25        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     26
    2127        /*if 2d: Error*/
    22         if (iomodel->dim==2) _error_("2d meshes not supported yet");
     28        if (dim==2) _error_("2d meshes not supported yet");
    2329
    2430        /*Recover pointer: */
    … …  
    3339        CreateSingleNodeToElementConnectivity(iomodel);
    3440
    35         for (i=0;i<iomodel->numberofvertices;i++){
     41        for (i=0;i<numberofvertices;i++){
    3642               
    3743                if((iomodel->my_vertices[i]==1)){

  • issm/trunk/src/c/modules/ModelProcessorx/Melting/CreateNodesMelting.cpp

    r9340 r9356  
    1818        int i;
    1919        bool continuous_galerkin=true;
     20        int    numberofvertices;
    2021
    2122        /*DataSets: */
    2223        Nodes*    nodes = NULL;
     24       
     25        /*Fetch parameters: */
     26        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    2327
    2428        /*Recover pointer: */
    … …  
    3943        iomodel->FetchData(&iomodel->nodeonwater,NULL,NULL,NodeOnWaterEnum);
    4044
    41         for (i=0;i<iomodel->numberofvertices;i++){
     45        for (i=0;i<numberofvertices;i++){
    4246
    4347                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/Melting/UpdateElementsMelting.cpp

    r9343 r9356  
    1616void    UpdateElementsMelting(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    dim;
     19        int    numberofelements;
     20
     21        /*Fetch parameters: */
     22        iomodel->parameters->FindParam(&dim,DimEnum);
     23        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     24
    1825        /*Now, is the model 3d? otherwise, do nothing: */
    19         if (iomodel->dim==2)return;
     26        if (dim==2)return;
    2027
    2128        /*Fetch data needed: */
    … …  
    2431        /*Update elements: */
    2532        int counter=0;
    26         for(int i=0;i<iomodel->numberofelements;i++){
     33        for(int i=0;i<numberofelements;i++){
    2734                if(iomodel->my_elements[i]){
    2835                        Element* element=(Element*)elements->GetObjectByOffset(counter);

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

    r9340 r9356  
    2020        int i;
    2121        int analysis_type;
     22        int dim;
     23        int verbose;
    2224       
    2325        /*output: */
    … …  
    3032        Parameters  *parameters  = NULL;
    3133
     34       
    3235        /*Initialize IoModel from input file*/
    3336        IoModel* iomodel = new IoModel(IOMODEL);
    34         SetVerbosityLevel(iomodel->verbose);
     37
     38        /*Fetch parameters: */
     39        iomodel->parameters->FindParam(&dim,DimEnum);
     40        iomodel->parameters->FindParam(&verbose,VerboseEnum);
     41       
     42        SetVerbosityLevel(verbose);
    3543
    3644        for(i=0;i<nummodels;i++){
    … …  
    3947
    4048                /*Hack for trasient runs (to be improved)*/
    41                 if(solution_type==TransientSolutionEnum && analysis_type==ThermalAnalysisEnum && iomodel->dim==2) continue;
    42                 if(solution_type==TransientSolutionEnum && analysis_type==MeltingAnalysisEnum && iomodel->dim==2) continue;
     49                if(solution_type==TransientSolutionEnum && analysis_type==ThermalAnalysisEnum && dim==2) continue;
     50                if(solution_type==TransientSolutionEnum && analysis_type==MeltingAnalysisEnum && dim==2) continue;
    4351       
    4452                _printf_(VerboseMProcessor(),"   create datasets for analysis %s\n",EnumToStringx(analysis_type));

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

    r9340 r9356  
    3636
    3737        /*as many nodes as there are vertices */
     38        int    numberofvertices;
    3839
    3940        /*output: */
    4041        bool* my_nodes=NULL;
    4142
    42         my_nodes=(bool*)xmalloc(iomodel->numberofvertices*sizeof(bool));
    43         memcpy(my_nodes,my_vertices,iomodel->numberofvertices*sizeof(bool));
     43        /*Fetch parameters: */
     44        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     45
     46        my_nodes=(bool*)xmalloc(numberofvertices*sizeof(bool));
     47        memcpy(my_nodes,my_vertices,numberofvertices*sizeof(bool));
    4448
    4549        /*Assign output pointers:*/
    … …  
    5054void  DiscontinuousGalerkinNodesPartitioning(bool** pmy_nodes,bool* my_elements, int* my_vertices, IoModel* iomodel){
    5155
     56        int    numberofelements;
     57
     58        /*Fetch parameters: */
     59        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     60
    5261        /*each element has it own nodes (as many as vertices) + additional nodes from neighbouring elements for each edge. This yields to a very different partition for
    5362         * the nodes and the vertices. The vertices are similar to continuous galerkin, but the nodes partitioning involves edges, which mess up sorting of
    … …  
    5564       
    5665        int i,j;
     66        int    dim;
    5767
    5868        /*output: */
    … …  
    6373        double  e1,e2;
    6474        int     pos;
     75        int     numberofedges;
     76
     77        /*Fetch parameters: */
     78        iomodel->parameters->FindParam(&dim,DimEnum);
    6579
    6680        /*Build discontinuous node partitioning
    … …  
    7488
    7589        /*Allocate*/
    76         my_nodes=(bool*)xcalloc(3*iomodel->numberofelements,sizeof(int));
     90        my_nodes=(bool*)xcalloc(3*numberofelements,sizeof(int));
    7791
    7892        /*First: add all the nodes of all the elements belonging to this cpu*/
    79         if (iomodel->dim==2){
    80                 for (i=0;i<iomodel->numberofelements;i++){
     93        if (dim==2){
     94                for (i=0;i<numberofelements;i++){
    8195                        if (my_elements[i]){
    8296                                my_nodes[3*i+0]=1;
    … …  
    93107
    94108        /*Get edges and elements*/
    95         iomodel->FetchData(&iomodel->edges,&iomodel->numberofedges,&cols,EdgesEnum);
     109        iomodel->FetchData(&iomodel->edges,&numberofedges,&cols,EdgesEnum);
    96110        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
    97111        if (cols!=4) _error_("field edges should have 4 columns");
    98112
    99113        /*!All elements have been partitioned above, only create elements for this CPU: */
    100         for (i=0;i<iomodel->numberofedges;i++){
     114        for (i=0;i<numberofedges;i++){
    101115
    102116                /*Get left and right elements*/

  • issm/trunk/src/c/modules/ModelProcessorx/Prognostic/CreateConstraintsPrognostic.cpp

    r9340 r9356  
    1010void    CreateConstraintsPrognostic(Constraints** pconstraints, IoModel* iomodel){
    1111
     12        int prognostic_DG;
     13       
     14        /*Fetch parameters: */
     15        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
     16
    1217        /*Output*/
    1318        Constraints *constraints = NULL;
    … …  
    2025
    2126        /*Do not add constraints in DG, they are weakly imposed*/
    22         if(!iomodel->prognostic_DG){
     27        if(!prognostic_DG){
    2328                IoModelToConstraintsx(constraints,iomodel,SpcthicknessEnum,PrognosticAnalysisEnum);
    2429        }

  • issm/trunk/src/c/modules/ModelProcessorx/Prognostic/CreateLoadsPrognostic.cpp

    r9340 r9356  
    1818        int penpair_ids[2];
    1919        int count=0;
     20        int prognostic_DG;
     21        int numberofedges;
     22        int numpenalties;
    2023
    2124        /*DataSet*/
    2225        Loads*    loads    = NULL;
     26       
     27        /*Fetch parameters: */
     28        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
    2329
    2430        /*Recover pointer: */
    … …  
    2935
    3036        /*Loads only in DG*/
    31         if (iomodel->prognostic_DG){
     37        if (prognostic_DG){
    3238
    3339                /*Get edges and elements*/
    34                 iomodel->FetchData(&iomodel->edges,&iomodel->numberofedges,NULL,EdgesEnum);
     40                iomodel->FetchData(&iomodel->edges,&numberofedges,NULL,EdgesEnum);
    3541                iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
    3642                iomodel->FetchData(&iomodel->thickness,NULL,NULL,ThicknessEnum);
    3743
    3844                /*First load data:*/
    39                 for (i=0;i<iomodel->numberofedges;i++){
     45                for (i=0;i<numberofedges;i++){
    4046
    4147                        /*Get left and right elements*/
    … …  
    5763
    5864        /*Create Penpair for penalties: */
    59         iomodel->FetchData(&iomodel->penalties,&iomodel->numpenalties,NULL,PenaltiesEnum);
     65        iomodel->FetchData(&iomodel->penalties,&numpenalties,NULL,PenaltiesEnum);
    6066        iomodel->FetchData(&iomodel->nodeonbed,NULL,NULL,NodeOnBedEnum);
    6167
    62         for(i=0;i<iomodel->numpenalties;i++){
     68        for(i=0;i<numpenalties;i++){
    6369
    6470                if(iomodel->my_vertices[(int)iomodel->penalties[2*i+0]-1]){

  • issm/trunk/src/c/modules/ModelProcessorx/Prognostic/CreateNodesPrognostic.cpp

    r9340 r9356  
    2121        int  io_index;
    2222        bool continuous_galerkin=true;
     23        int  dim;
     24        int    numberofelements;
     25        int    numberofvertices;
     26        int    prognostic_DG;
    2327
    2428        /*DataSets: */
    2529        Nodes*    nodes = NULL;
     30
     31        /*Fetch parameters: */
     32        iomodel->parameters->FindParam(&dim,DimEnum);
     33        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     34        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     35        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
    2636
    2737        /*Recover pointer: */
    … …  
    3242
    3343        /*Create partition of nodes: */
    34         if(iomodel->prognostic_DG) continuous_galerkin=false;
     44        if(prognostic_DG) continuous_galerkin=false;
    3545        NodesPartitioning(&iomodel->my_nodes,iomodel->my_elements,iomodel->my_vertices,iomodel,continuous_galerkin);
    3646
    3747        /*Check in 3d*/
    38         if(iomodel->prognostic_DG && iomodel->dim==3) _error_("DG 3d not implemented yet");
     48        if(prognostic_DG && dim==3) _error_("DG 3d not implemented yet");
    3949
    4050        /*First fetch data: */
    … …  
    5060
    5161                /*Build Nodes dataset (Continuous Galerkin)*/
    52                 for (i=0;i<iomodel->numberofvertices;i++){
     62                for (i=0;i<numberofvertices;i++){
    5363
    5464                        if(iomodel->my_vertices[i]){
    … …  
    6373
    6474                /*Build Nodes dataset -> 3 for each element (Discontinuous Galerkin)*/
    65                 for (i=0;i<iomodel->numberofelements;i++){
     75                for (i=0;i<numberofelements;i++){
    6676                        for (j=0;j<3;j++){
    6777
    … …  
    7181                                        vertex_id=(int)*(iomodel->elements+3*i+j); //(Matlab indexing)
    7282                                        io_index=vertex_id-1;                      //(C indexing)
    73                                         _assert_(vertex_id>0 && vertex_id<=iomodel->numberofvertices);
     83                                        _assert_(vertex_id>0 && vertex_id<=numberofvertices);
    7484
    7585                                        //Compute Node id

  • issm/trunk/src/c/modules/ModelProcessorx/Prognostic/UpdateElementsPrognostic.cpp

    r9343 r9356  
    1616void    UpdateElementsPrognostic(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    dim;
     19        int    numberofelements;
     20        int    prognostic_DG;
     21
    1822        /*Fetch data needed: */
     23        iomodel->parameters->FindParam(&dim,DimEnum);
    1924        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
     25        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     26        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
    2027
    2128        /*Update elements: */
    2229        int counter=0;
    23         for(int i=0;i<iomodel->numberofelements;i++){
     30        for(int i=0;i<numberofelements;i++){
    2431                if(iomodel->my_elements[i]){
    2532                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    4249        iomodel->FetchDataToInput(elements,VyEnum);
    4350
    44         if(iomodel->prognostic_DG){
     51        if(prognostic_DG){
    4552                iomodel->FetchDataToInput(elements,SpcthicknessEnum); //for DG, we need the spc in the element
    4653        }
    4754       
    48         if (iomodel->dim==3){
     55        if (dim==3){
    4956                iomodel->FetchDataToInput(elements,ElementOnBedEnum);
    5057                iomodel->FetchDataToInput(elements,ElementOnSurfaceEnum);

  • issm/trunk/src/c/modules/ModelProcessorx/Qmu/CreateParametersQmu.cpp

    r9340 r9356  
    5252        int      m,n;
    5353        int      count;
     54        bool     qmu_analysis=false;
     55        char*    name=NULL;
     56        int      numberofresponses;
     57        int      numberofvertices;
    5458
    5559        /*}}}*/
    … …  
    5862        parameters=*pparameters;
    5963
    60         parameters->AddObject(new   BoolParam(QmuAnalysisEnum,iomodel->qmu_analysis));
    61         if(iomodel->qmu_analysis){
     64        /*recover parameters: */
     65        parameters->FindParam(&qmu_analysis,QmuAnalysisEnum);
     66        parameters->FindParam(&name,NameEnum);
     67        parameters->FindParam(&numberofresponses,NumberOfResponsesEnum);
     68        parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     69        if(qmu_analysis){
    6270
    6371                /*name of qmu input, error and output files:{{{1*/
    64                 qmuinname=(char*)xmalloc((strlen(iomodel->name)+strlen(".qmu.in")+1)*sizeof(char));
    65                 sprintf(qmuinname,"%s%s",iomodel->name,".qmu.in");
     72                qmuinname=(char*)xmalloc((strlen(name)+strlen(".qmu.in")+1)*sizeof(char));
     73                sprintf(qmuinname,"%s%s",name,".qmu.in");
    6674                parameters->AddObject(new   StringParam(QmuInNameEnum,qmuinname));
    6775
    68                 qmuoutname=(char*)xmalloc((strlen(iomodel->name)+strlen(".qmu.out")+1)*sizeof(char));
    69                 sprintf(qmuoutname,"%s%s",iomodel->name,".qmu.out");
     76                qmuoutname=(char*)xmalloc((strlen(name)+strlen(".qmu.out")+1)*sizeof(char));
     77                sprintf(qmuoutname,"%s%s",name,".qmu.out");
    7078                parameters->AddObject(new   StringParam(QmuOutNameEnum,qmuoutname));
    7179
    72                 qmuerrname=(char*)xmalloc((strlen(iomodel->name)+strlen(".qmu.err")+1)*sizeof(char));
    73                 sprintf(qmuerrname,"%s%s",iomodel->name,".qmu.err");
     80                qmuerrname=(char*)xmalloc((strlen(name)+strlen(".qmu.err")+1)*sizeof(char));
     81                sprintf(qmuerrname,"%s%s",name,".qmu.err");
    7482                parameters->AddObject(new   StringParam(QmuErrNameEnum,qmuerrname));
    7583                /*}}}*/
    … …  
    8694                /*Ok, we have all the response descriptors. Build a parameter with it: */
    8795                parameters->AddObject(new StringArrayParam(ResponsedescriptorsEnum,responsedescriptors,numresponsedescriptors));
    88                 parameters->AddObject(new    IntParam(QmuNumberOfResponsesEnum,iomodel->numberofresponses));
     96                parameters->AddObject(new    IntParam(QmuNumberOfResponsesEnum,numberofresponses));
    8997                /*}}}*/
    9098                /*Deal with partitioning: {{{1*/
    9199                /*partition vertices in iomodel->qmu_npart parts, unless a partition is already present: */
    92                 parameters->AddObject(new    IntParam(QmuNPartEnum,iomodel->qmu_npart));
    93100               
    94101                iomodel->FetchData(&dpart,NULL,NULL,PartEnum);
    … …  
    99106                        ElementsAndVerticesPartitioning(&iomodel->my_elements,&iomodel->my_vertices,iomodel);
    100107
    101                         dpart=(double*)xmalloc(iomodel->numberofvertices*sizeof(double));
    102                         for(i=0;i<iomodel->numberofvertices;i++)dpart[i]=iomodel->my_vertices[i];
    103                 }
    104                 parameters->AddObject(new DoubleVecParam(QmuPartEnum,dpart,iomodel->numberofvertices));
     108                        dpart=(double*)xmalloc(numberofvertices*sizeof(double));
     109                        for(i=0;i<numberofvertices;i++)dpart[i]=iomodel->my_vertices[i];
     110                }
     111                parameters->AddObject(new DoubleVecParam(QmuPartEnum,dpart,numberofvertices));
    105112                /*}}}*/
    106113                /*Deal with data needed because of qmu variables: {{{1*/
    … …  
    116123
    117124                                /*Convert units: */
    118                                 UnitConversion(dakota_parameter,iomodel->numberofvertices,ExtToIuEnum,StringToEnumx(tag));
     125                                UnitConversion(dakota_parameter,numberofvertices,ExtToIuEnum,StringToEnumx(tag));
    119126
    120127                                /*Add to parameters: */
    121                                 parameters->AddObject(new DoubleVecParam(StringToEnumx(tag),dakota_parameter,iomodel->numberofvertices));
     128                                parameters->AddObject(new DoubleVecParam(StringToEnumx(tag),dakota_parameter,numberofvertices));
    122129                               
    123130                                /*Free ressources:*/
    … …  
    187194                }
    188195                /*}}}*/
    189                 parameters->AddObject(new   BoolParam(QmuSaveFemmodelEnum,iomodel->qmu_save_femmodel));
    190196                /*Free data: {{{1*/
    191197                for(i=0;i<numresponsedescriptors;i++){
    … …  
    208214                xfree((void**)&qmuerrname);
    209215                xfree((void**)&qmuoutname);
    210                 /*}}}*/
    211         } //if(iomodel->qmu_analysis)
     216                xfree((void**)&name);
     217                /*}}}*/
     218        } //if(qmu_analysis)
    212219
    213220        /*Assign output pointer: */

  • issm/trunk/src/c/modules/ModelProcessorx/SurfaceSlope/CreateNodesSurfaceSlope.cpp

    r9340 r9356  
    1818        int i;
    1919        bool continuous_galerkin=true;
     20        int    numberofvertices;
    2021
    2122        /*DataSets: */
    2223        Nodes*    nodes = NULL;
     24
     25        /*Fetch parameters: */
     26        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
    2327
    2428        /*Recover pointer: */
    … …  
    3943        iomodel->FetchData(&iomodel->nodeonwater,NULL,NULL,NodeOnWaterEnum);
    4044
    41         for (i=0;i<iomodel->numberofvertices;i++){
     45        for (i=0;i<numberofvertices;i++){
    4246
    4347                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/SurfaceSlope/UpdateElementsSurfaceSlope.cpp

    r9343 r9356  
    1616void    UpdateElementsSurfaceSlope(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    dim;
     19        int    numberofelements;
     20
    1821        /*Fetch data needed: */
     22        iomodel->parameters->FindParam(&dim,DimEnum);
    1923        iomodel->FetchData(&iomodel->elements,NULL,NULL,ElementsEnum);
     24        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
    2025
    2126        /*Update elements: */
    2227        int counter=0;
    23         for(int i=0;i<iomodel->numberofelements;i++){
     28        for(int i=0;i<numberofelements;i++){
    2429                if(iomodel->my_elements[i]){
    2530                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    3338        iomodel->FetchDataToInput(elements,ElementOnWaterEnum);
    3439       
    35         if (iomodel->dim==3){
     40        if (dim==3){
    3641                iomodel->FetchDataToInput(elements,ElementOnBedEnum);
    3742                iomodel->FetchDataToInput(elements,ElementOnSurfaceEnum);

  • issm/trunk/src/c/modules/ModelProcessorx/Thermal/CreateConstraintsThermal.cpp

    r9340 r9356  
    1717        int i;
    1818        int count;
     19        int    dim;
    1920       
    2021        /*Output*/
    2122        Constraints* constraints = NULL;
     23
     24        /*Fetch parameters: */
     25        iomodel->parameters->FindParam(&dim,DimEnum);
    2226
    2327        /*Recover pointer: */
    … …  
    2832
    2933        /*Only 3d mesh supported*/
    30         if (iomodel->dim==3){
     34        if (dim==3){
    3135                IoModelToConstraintsx(constraints,iomodel,SpctemperatureEnum,ThermalAnalysisEnum);
    3236        }

  • issm/trunk/src/c/modules/ModelProcessorx/Thermal/CreateLoadsThermal.cpp

    r9340 r9356  
    1515        /*Intermediary*/
    1616        int i;
     17        int    dim;
     18        int    numberofvertices;
    1719
    1820        /*DataSet*/
    … …  
    2325        loads=*ploads;
    2426
     27        /*Fetch parameters: */
     28        iomodel->parameters->FindParam(&dim,DimEnum);
     29        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     30
    2531        /*Create loads if they do not exist yet*/
    2632        if(!loads) loads = new Loads(LoadsEnum);
    2733
    2834        /*return if 2d mesh*/
    29         if (iomodel->dim==2) _error_("2d meshes not supported yet");
     35        if (dim==2) _error_("2d meshes not supported yet");
    3036
    3137        //create penalties for nodes: no node can have a temperature over the melting point
    … …  
    3440        CreateSingleNodeToElementConnectivity(iomodel);
    3541
    36         for (i=0;i<iomodel->numberofvertices;i++){
     42        for (i=0;i<numberofvertices;i++){
    3743       
    3844                /*keep only this partition's nodes:*/

  • issm/trunk/src/c/modules/ModelProcessorx/Thermal/CreateNodesThermal.cpp

    r9340 r9356  
    1818        int i;
    1919        bool continuous_galerkin=true;
     20        int    numberofvertices;
    2021
    2122        /*DataSets: */
    2223        Nodes*    nodes = NULL;
    23        
     24
     25        /*Fetch parameters: */
     26        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     27
    2428        /*Recover pointer: */
    2529        nodes=*pnodes;
    … …  
    3943        iomodel->FetchData(&iomodel->nodeonwater,NULL,NULL,NodeOnWaterEnum);
    4044
    41         for (i=0;i<iomodel->numberofvertices;i++){
     45        for (i=0;i<numberofvertices;i++){
    4246
    4347                if(iomodel->my_vertices[i]){

  • issm/trunk/src/c/modules/ModelProcessorx/Thermal/UpdateElementsThermal.cpp

    r9343 r9356  
    1616void    UpdateElementsThermal(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
     18        int    dim;
     19        int    numberofelements;
     20        bool   qmu_analysis;
     21
     22        /*Fetch parameters: */
     23        iomodel->parameters->FindParam(&dim,DimEnum);
     24        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     25        iomodel->parameters->FindParam(&qmu_analysis,QmuAnalysisEnum);
     26
    1827        /*Now, is the model 3d? otherwise, do nothing: */
    19         if (iomodel->dim==2)return;
     28        if (dim==2)return;
    2029
    2130        /*Fetch data needed: */
    … …  
    2433        /*Update elements: */
    2534        int counter=0;
    26         for(int i=0;i<iomodel->numberofelements;i++){
     35        for(int i=0;i<numberofelements;i++){
    2736                if(iomodel->my_elements[i]){
    2837                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    5362        iomodel->FetchDataToInput(elements,VzEnum);
    5463       
    55         if(iomodel->qmu_analysis)elements->InputDuplicate(TemperatureEnum,QmuTemperatureEnum);
     64        if(qmu_analysis)elements->InputDuplicate(TemperatureEnum,QmuTemperatureEnum);
    5665
    5766        /*Free data: */

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

    r9291 r9356  
    143143        else if (strcmp(name,"Adjointz")==0) return AdjointzEnum;
    144144        else if (strcmp(name,"Adjointp")==0) return AdjointpEnum;
    145         else if (strcmp(name,"ArtDiff")==0) return ArtDiffEnum;
     145        else if (strcmp(name,"Fake35")==0) return Fake35Enum;
    146146        else if (strcmp(name,"Bed")==0) return BedEnum;
    147147        else if (strcmp(name,"Bathymetry")==0) return BathymetryEnum;
    … …  
    153153        else if (strcmp(name,"Fake2")==0) return Fake2Enum;
    154154        else if (strcmp(name,"Constant")==0) return ConstantEnum;
    155         else if (strcmp(name,"NumControls")==0) return NumControlsEnum;
     155        else if (strcmp(name,"Fake37")==0) return Fake37Enum;
    156156        else if (strcmp(name,"ControlType")==0) return ControlTypeEnum;
    157157        else if (strcmp(name,"Converged")==0) return ConvergedEnum;
    … …  
    348348        else if (strcmp(name,"QmuInName")==0) return QmuInNameEnum;
    349349        else if (strcmp(name,"QmuMassFluxSegments")==0) return QmuMassFluxSegmentsEnum;
    350         else if (strcmp(name,"QmuNPart")==0) return QmuNPartEnum;
     350        else if (strcmp(name,"Fake39")==0) return Fake39Enum;
    351351        else if (strcmp(name,"QmuOutName")==0) return QmuOutNameEnum;
    352352        else if (strcmp(name,"QmuPart")==0) return QmuPartEnum;
    … …  
    370370        else if (strcmp(name,"Fset")==0) return FsetEnum;
    371371        else if (strcmp(name,"Sset")==0) return SsetEnum;
    372         else if (strcmp(name,"GroundingLineMigration")==0) return GroundingLineMigrationEnum;
     372        else if (strcmp(name,"GroundinglineMigration")==0) return GroundinglineMigrationEnum;
    373373        else if (strcmp(name,"Yts")==0) return YtsEnum;
    374374        else if (strcmp(name,"TriangleInterp")==0) return TriangleInterpEnum;
    … …  
    399399        else if (strcmp(name,"VelAbsGradient")==0) return VelAbsGradientEnum;
    400400        else if (strcmp(name,"DatasetInput")==0) return DatasetInputEnum;
    401         else if (strcmp(name,"NumResponses")==0) return NumResponsesEnum;
     401        else if (strcmp(name,"Fake38")==0) return Fake38Enum;
    402402        else if (strcmp(name,"StepResponses")==0) return StepResponsesEnum;
    403403        else if (strcmp(name,"IntMatParam")==0) return IntMatParamEnum;
    … …  
    429429        else if (strcmp(name,"NodeOnHutter")==0) return NodeOnHutterEnum;
    430430        else if (strcmp(name,"Z")==0) return ZEnum;
    431         else if (strcmp(name,"GlMigration")==0) return GlMigrationEnum;
     431        else if (strcmp(name,"Fake36")==0) return Fake36Enum;
    432432        else if (strcmp(name,"Riftinfo")==0) return RiftinfoEnum;
    433433        else if (strcmp(name,"ElementOnIceSheet")==0) return ElementOnIceSheetEnum;

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

    r9320 r9356  
    4949Penta::Penta(int penta_id, int penta_sid, int index, IoModel* iomodel,int nummodels)
    5050        :PentaRef(nummodels)
    51         ,PentaHook(nummodels,index+1,iomodel->numberofelements+1) //index+1: matice id, iomodel->numberofelements+1: matpar id
     51        ,PentaHook(nummodels,index+1,iomodel) //index+1: matice id, iomodel->numberofelements+1: matpar id
    5252                                                                      { //i is the element index
    5353
    … …  
    18811881        thermalconductivity=matpar->GetThermalConductivity();
    18821882        this->parameters->FindParam(&dt,DtEnum);
    1883         this->parameters->FindParam(&artdiff,ArtDiffEnum);
     1883        this->parameters->FindParam(&artdiff,ArtificialDiffusivityEnum);
    18841884        this->parameters->FindParam(&epsvel,EpsvelEnum);
    18851885        Input* pressure_input=inputs->GetInput(PressureEnum);      _assert_(pressure_input);
    … …  
    21552155        thermalconductivity=matpar->GetThermalConductivity();
    21562156        this->parameters->FindParam(&dt,DtEnum);
    2157         this->parameters->FindParam(&artdiff,ArtDiffEnum);
     2157        this->parameters->FindParam(&artdiff,ArtificialDiffusivityEnum);
    21582158        this->parameters->FindParam(&epsvel,EpsvelEnum);
    21592159        Input* vx_input=inputs->GetInput(VxEnum);      _assert_(vx_input);
    … …  
    33183318        thermalconductivity=matpar->GetThermalConductivity();
    33193319        this->parameters->FindParam(&dt,DtEnum);
    3320         this->parameters->FindParam(&artdiff,ArtDiffEnum);
     3320        this->parameters->FindParam(&artdiff,ArtificialDiffusivityEnum);
    33213321        Input* vx_input=inputs->GetInput(VxEnum); _assert_(vx_input);
    33223322        Input* vy_input=inputs->GetInput(VyEnum); _assert_(vy_input);
    … …  
    35873587        thermalconductivity=matpar->GetThermalConductivity();
    35883588        this->parameters->FindParam(&dt,DtEnum);
    3589         this->parameters->FindParam(&artdiff,ArtDiffEnum);
     3589        this->parameters->FindParam(&artdiff,ArtificialDiffusivityEnum);
    35903590        Input* vx_input=inputs->GetInput(VxEnum); _assert_(vx_input);
    35913591        Input* vy_input=inputs->GetInput(VyEnum); _assert_(vy_input);
    … …  
    44674467        int        *responses = NULL;
    44684468        int         num_responses,resp;
    4469         this->parameters->FindParam(&num_responses,NumResponsesEnum);
     4469        this->parameters->FindParam(&num_responses,NumCmResponsesEnum);
    44704470        this->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum);
    44714471
    … …  
    47954795
    47964796        /*retrieve some parameters: */
    4797         this->parameters->FindParam(&num_controls,NumControlsEnum);
     4797        this->parameters->FindParam(&num_controls,NumControlTypeEnum);
    47984798        this->parameters->FindParam(&control_type,NULL,ControlTypeEnum);
    47994799
    … …  
    48534853        double time;
    48544854        TransientInput* transientinput=NULL;
    4855                
     4855
     4856        int    numberofvertices;
     4857        int    numberofelements;
     4858        double yts;
     4859
     4860        /*Fetch parameters: */
     4861        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     4862        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     4863        iomodel->parameters->FindParam(&yts,YtsEnum);
     4864
    48564865        /*Branch on type of vector: nodal or elementary: */
    48574866        if(vector_type==1){ //nodal vector
    … …  
    48634872
    48644873                /*Are we in transient or static? */
    4865                 if(M==iomodel->numberofvertices){
     4874                if(M==numberofvertices){
    48664875
    48674876                        /*create input values: */
    … …  
    48744883                        this->inputs->AddInput(new PentaVertexInput(vector_enum,nodeinputs));
    48754884                }
    4876                 else if(M==iomodel->numberofvertices+1){
     4885                else if(M==numberofvertices+1){
    48774886                        /*create transient input: */
    48784887                        for(t=0;t<N;t++){ //N is the number of times
    … …  
    48884897
    48894898                                /*time? :*/
    4890                                 time=(double)vector[(M-1)*N+t]*iomodel->yts;
     4899                                time=(double)vector[(M-1)*N+t]*yts;
    48914900
    48924901                                if(t==0)transientinput=new TransientInput(vector_enum);
    … …  
    48954904                        this->inputs->AddInput(transientinput);
    48964905                }
    4897                 else _error_("nodal vector is either numberofnodes (%i), or numberofnodes+1 long. Field provided is %i long. Enum %s",iomodel->numberofvertices,M,EnumToStringx(vector_enum));
     4906                else _error_("nodal vector is either numberofnodes (%i), or numberofnodes+1 long. Field provided is %i long. Enum %s",numberofvertices,M,EnumToStringx(vector_enum));
    48984907        }
    48994908        else if(vector_type==2){ //element vector
    49004909                /*Are we in transient or static? */
    4901                 if(M==iomodel->numberofelements){
     4910                if(M==numberofelements){
    49024911
    49034912                        /*static mode: create an input out of the element value: */
    … …  
    51535162        double  cmmaxinputs[6];
    51545163
     5164        double  yts;
     5165        int     control_analysis;
     5166        int     num_control_type;
     5167        int     num_cm_responses;
     5168
     5169        /*Fetch parameters: */
     5170        iomodel->parameters->FindParam(&yts,YtsEnum);
     5171        iomodel->parameters->FindParam(&control_analysis,ControlAnalysisEnum);
     5172        iomodel->parameters->FindParam(&num_control_type,NumControlTypeEnum);
     5173        iomodel->parameters->FindParam(&num_cm_responses,NumCmResponsesEnum);
     5174
    51555175        /*Checks if debuging*/
    51565176        /*{{{2*/
    … …  
    51655185       
    51665186        /*Control Inputs*/
    5167         if (iomodel->control_analysis && iomodel->control_type){
    5168                 for(i=0;i<iomodel->num_control_type;i++){
     5187        if (control_analysis && iomodel->control_type){
     5188                for(i=0;i<num_control_type;i++){
    51695189                        switch((int)iomodel->control_type[i]){
    51705190                                case DhdtEnum:
    51715191                                        if (iomodel->dhdt){
    5172                                                 for(j=0;j<6;j++)nodeinputs[j]=iomodel->dhdt[penta_vertex_ids[j]-1]/iomodel->yts;
    5173                                                 for(j=0;j<6;j++)cmmininputs[j]=iomodel->cm_min[(penta_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
    5174                                                 for(j=0;j<6;j++)cmmaxinputs[j]=iomodel->cm_max[(penta_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
     5192                                                for(j=0;j<6;j++)nodeinputs[j]=iomodel->dhdt[penta_vertex_ids[j]-1]/yts;
     5193                                                for(j=0;j<6;j++)cmmininputs[j]=iomodel->cm_min[(penta_vertex_ids[j]-1)*num_control_type+i]/yts;
     5194                                                for(j=0;j<6;j++)cmmaxinputs[j]=iomodel->cm_max[(penta_vertex_ids[j]-1)*num_control_type+i]/yts;
    51755195                                                this->inputs->AddInput(new ControlInput(DhdtEnum,PentaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    51765196                                        }
    … …  
    51785198                                case VxEnum:
    51795199                                        if (iomodel->vx){
    5180                                                 for(j=0;j<6;j++)nodeinputs[j]=iomodel->vx[penta_vertex_ids[j]-1]/iomodel->yts;
    5181                                                 for(j=0;j<6;j++)cmmininputs[j]=iomodel->cm_min[(penta_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
    5182                                                 for(j=0;j<6;j++)cmmaxinputs[j]=iomodel->cm_max[(penta_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
     5200                                                for(j=0;j<6;j++)nodeinputs[j]=iomodel->vx[penta_vertex_ids[j]-1]/yts;
     5201                                                for(j=0;j<6;j++)cmmininputs[j]=iomodel->cm_min[(penta_vertex_ids[j]-1)*num_control_type+i]/yts;
     5202                                                for(j=0;j<6;j++)cmmaxinputs[j]=iomodel->cm_max[(penta_vertex_ids[j]-1)*num_control_type+i]/yts;
    51835203                                                this->inputs->AddInput(new ControlInput(VxEnum,PentaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    51845204                                        }
    … …  
    51865206                                case VyEnum:
    51875207                                        if (iomodel->vy){
    5188                                                 for(j=0;j<6;j++)nodeinputs[j]=iomodel->vy[penta_vertex_ids[j]-1]/iomodel->yts;
    5189                                                 for(j=0;j<6;j++)cmmininputs[j]=iomodel->cm_min[(penta_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
    5190                                                 for(j=0;j<6;j++)cmmaxinputs[j]=iomodel->cm_max[(penta_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
     5208                                                for(j=0;j<6;j++)nodeinputs[j]=iomodel->vy[penta_vertex_ids[j]-1]/yts;
     5209                                                for(j=0;j<6;j++)cmmininputs[j]=iomodel->cm_min[(penta_vertex_ids[j]-1)*num_control_type+i]/yts;
     5210                                                for(j=0;j<6;j++)cmmaxinputs[j]=iomodel->cm_max[(penta_vertex_ids[j]-1)*num_control_type+i]/yts;
    51915211                                                this->inputs->AddInput(new ControlInput(VyEnum,PentaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    51925212                                        }
    … …  
    51955215                                        if (iomodel->drag_coefficient){
    51965216                                                for(j=0;j<6;j++)nodeinputs[j]=iomodel->drag_coefficient[penta_vertex_ids[j]-1];
    5197                                                 for(j=0;j<6;j++)cmmininputs[j]=iomodel->cm_min[(penta_vertex_ids[j]-1)*iomodel->num_control_type+i];
    5198                                                 for(j=0;j<6;j++)cmmaxinputs[j]=iomodel->cm_max[(penta_vertex_ids[j]-1)*iomodel->num_control_type+i];
     5217                                                for(j=0;j<6;j++)cmmininputs[j]=iomodel->cm_min[(penta_vertex_ids[j]-1)*num_control_type+i];
     5218                                                for(j=0;j<6;j++)cmmaxinputs[j]=iomodel->cm_max[(penta_vertex_ids[j]-1)*num_control_type+i];
    51995219                                                this->inputs->AddInput(new ControlInput(DragCoefficientEnum,PentaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    52005220                                        }
    … …  
    52445264                /*Create inputs and add to DataSetInput*/
    52455265                DatasetInput* datasetinput=new DatasetInput(WeightsEnum);
    5246                 for(i=0;i<iomodel->num_cm_responses;i++){
    5247                         for(j=0;j<6;j++)nodeinputs[j]=iomodel->weights[(penta_vertex_ids[j]-1)*iomodel->num_cm_responses+i];
     5266                for(i=0;i<num_cm_responses;i++){
     5267                        for(j=0;j<6;j++)nodeinputs[j]=iomodel->weights[(penta_vertex_ids[j]-1)*num_cm_responses+i];
    52485268                        datasetinput->inputs->AddObject(new PentaVertexInput(WeightsEnum,nodeinputs));
    52495269                }
    … …  
    74937513        int     penta_vertex_ids[6];
    74947514        double  nodeinputs[6];
     7515        double  yts;
     7516        int     prognostic_DG;
     7517        bool    qmu_analysis;
     7518        int     isstokes;
     7519        double  beta,heatcapacity,referencetemperature,meltingpoint,latentheat;
     7520
     7521        /*Fetch parameters: */
     7522        iomodel->parameters->FindParam(&yts,YtsEnum);
     7523        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
     7524        iomodel->parameters->FindParam(&qmu_analysis,QmuAnalysisEnum);
     7525        iomodel->parameters->FindParam(&isstokes,IsstokesEnum);
     7526        iomodel->parameters->FindParam(&beta,BetaEnum);
     7527        iomodel->parameters->FindParam(&heatcapacity,HeatcapacityEnum);
     7528        iomodel->parameters->FindParam(&referencetemperature,ReferencetemperatureEnum);
     7529        iomodel->parameters->FindParam(&meltingpoint,MeltingpointEnum);
     7530        iomodel->parameters->FindParam(&latentheat,LatentheatEnum);
     7531
    74957532
    74967533        /*Checks if debuging*/
    … …  
    75007537
    75017538        /*Recover element type*/
    7502         if ((analysis_type==PrognosticAnalysisEnum || analysis_type==BalancethicknessAnalysisEnum) && iomodel->prognostic_DG){
     7539        if ((analysis_type==PrognosticAnalysisEnum || analysis_type==BalancethicknessAnalysisEnum) && prognostic_DG){
    75037540                /*P1 Discontinuous Galerkin*/
    75047541                penta_type=P1DGEnum;
    … …  
    75337570                        /*default vx,vy and vz: either observation or 0 */
    75347571                        if(!iomodel->vx){
    7535                                 if (iomodel->vx_obs) for(i=0;i<6;i++)nodeinputs[i]=iomodel->vx_obs[penta_vertex_ids[i]-1]/iomodel->yts;
     7572                                if (iomodel->vx_obs) for(i=0;i<6;i++)nodeinputs[i]=iomodel->vx_obs[penta_vertex_ids[i]-1]/yts;
    75367573                                else                 for(i=0;i<6;i++)nodeinputs[i]=0;
    75377574                                this->inputs->AddInput(new PentaVertexInput(VxEnum,nodeinputs));
    75387575                                this->inputs->AddInput(new PentaVertexInput(VxPicardEnum,nodeinputs));
    7539                                 if(iomodel->qmu_analysis) this->inputs->AddInput(new PentaVertexInput(QmuVxEnum,nodeinputs));
     7576                                if(qmu_analysis) this->inputs->AddInput(new PentaVertexInput(QmuVxEnum,nodeinputs));
    75407577                        }
    75417578                        if(!iomodel->vy){
    7542                                 if (iomodel->vy_obs) for(i=0;i<6;i++)nodeinputs[i]=iomodel->vy_obs[penta_vertex_ids[i]-1]/iomodel->yts;
     7579                                if (iomodel->vy_obs) for(i=0;i<6;i++)nodeinputs[i]=iomodel->vy_obs[penta_vertex_ids[i]-1]/yts;
    75437580                                else                 for(i=0;i<6;i++)nodeinputs[i]=0;
    75447581                                this->inputs->AddInput(new PentaVertexInput(VyEnum,nodeinputs));
    75457582                                this->inputs->AddInput(new PentaVertexInput(VyPicardEnum,nodeinputs));
    7546                                 if(iomodel->qmu_analysis) this->inputs->AddInput(new PentaVertexInput(QmuVyEnum,nodeinputs));
     7583                                if(qmu_analysis) this->inputs->AddInput(new PentaVertexInput(QmuVyEnum,nodeinputs));
    75477584                        }
    75487585                        if(!iomodel->vz){
    7549                                 if (iomodel->vz_obs) for(i=0;i<6;i++)nodeinputs[i]=iomodel->vz_obs[penta_vertex_ids[i]-1]/iomodel->yts;
     7586                                if (iomodel->vz_obs) for(i=0;i<6;i++)nodeinputs[i]=iomodel->vz_obs[penta_vertex_ids[i]-1]/yts;
    75507587                                else                 for(i=0;i<6;i++)nodeinputs[i]=0;
    75517588                                this->inputs->AddInput(new PentaVertexInput(VzEnum,nodeinputs));
    75527589                                this->inputs->AddInput(new PentaVertexInput(VzPicardEnum,nodeinputs));
    7553                                 if(iomodel->qmu_analysis) this->inputs->AddInput(new PentaVertexInput(QmuVzEnum,nodeinputs));
     7590                                if(qmu_analysis) this->inputs->AddInput(new PentaVertexInput(QmuVzEnum,nodeinputs));
    75547591                        }
    75557592                        if(!iomodel->pressure){
    75567593                                for(i=0;i<6;i++)nodeinputs[i]=0;
    7557                                 if(iomodel->qmu_analysis){
     7594                                if(qmu_analysis){
    75587595                                        this->inputs->AddInput(new PentaVertexInput(PressureEnum,nodeinputs));
    75597596                                        this->inputs->AddInput(new PentaVertexInput(QmuPressureEnum,nodeinputs));
    75607597                                }
    7561                                 if(iomodel->isstokes){
     7598                                if(isstokes){
    75627599                                        this->inputs->AddInput(new PentaVertexInput(PressureEnum,nodeinputs));
    75637600                                        this->inputs->AddInput(new PentaVertexInput(PressurePicardEnum,nodeinputs));
    … …  
    75677604                                /*Create VzPattyn and VzStokes Enums*/
    75687605                                if(iomodel->vz && iomodel->nodeonstokes){
    7569                                         for(i=0;i<6;i++) nodeinputs[i]=iomodel->vz[penta_vertex_ids[i]-1]/iomodel->yts*iomodel->nodeonstokes[penta_vertex_ids[i]-1];
     7606                                        for(i=0;i<6;i++) nodeinputs[i]=iomodel->vz[penta_vertex_ids[i]-1]/yts*iomodel->nodeonstokes[penta_vertex_ids[i]-1];
    75707607                                        this->inputs->AddInput(new PentaVertexInput(VzStokesEnum,nodeinputs));
    7571                                         for(i=0;i<6;i++) nodeinputs[i]=iomodel->vz[penta_vertex_ids[i]-1]/iomodel->yts*(1-iomodel->nodeonstokes[penta_vertex_ids[i]-1]);
     7608                                        for(i=0;i<6;i++) nodeinputs[i]=iomodel->vz[penta_vertex_ids[i]-1]/yts*(1-iomodel->nodeonstokes[penta_vertex_ids[i]-1]);
    75727609                                        this->inputs->AddInput(new PentaVertexInput(VzPattynEnum,nodeinputs));
    75737610                                }
    … …  
    75817618                                /*Create VzMacAyeal and VzStokes Enums*/
    75827619                                if(iomodel->vz && iomodel->nodeonstokes){
    7583                                         for(i=0;i<6;i++) nodeinputs[i]=iomodel->vz[penta_vertex_ids[i]-1]/iomodel->yts*iomodel->nodeonstokes[penta_vertex_ids[i]-1];
     7620                                        for(i=0;i<6;i++) nodeinputs[i]=iomodel->vz[penta_vertex_ids[i]-1]/yts*iomodel->nodeonstokes[penta_vertex_ids[i]-1];
    75847621                                        this->inputs->AddInput(new PentaVertexInput(VzStokesEnum,nodeinputs));
    7585                                         for(i=0;i<6;i++) nodeinputs[i]=iomodel->vz[penta_vertex_ids[i]-1]/iomodel->yts*(1-iomodel->nodeonstokes[penta_vertex_ids[i]-1]);
     7622                                        for(i=0;i<6;i++) nodeinputs[i]=iomodel->vz[penta_vertex_ids[i]-1]/yts*(1-iomodel->nodeonstokes[penta_vertex_ids[i]-1]);
    75867623                                        this->inputs->AddInput(new PentaVertexInput(VzMacAyealEnum,nodeinputs));
    75877624                                }
    … …  
    76107647                        if (iomodel->temperature && iomodel->waterfraction) {
    76117648                                for(i=0;i<6;i++){
    7612                                         if(iomodel->temperature[penta_vertex_ids[i]-1] < iomodel->meltingpoint-iomodel->beta*iomodel->pressure[penta_vertex_ids[i]-1]){
    7613                                                 nodeinputs[i]=iomodel->heatcapacity*(iomodel->temperature[penta_vertex_ids[i]-1]-iomodel->referencetemperature);
     7649                                        if(iomodel->temperature[penta_vertex_ids[i]-1] < meltingpoint-beta*iomodel->pressure[penta_vertex_ids[i]-1]){
     7650                                                nodeinputs[i]=heatcapacity*(iomodel->temperature[penta_vertex_ids[i]-1]-referencetemperature);
    76147651                                        }
    7615                                         else nodeinputs[i]=iomodel->heatcapacity*
    7616                                          (iomodel->meltingpoint-iomodel->beta*iomodel->pressure[penta_vertex_ids[i]-1]-iomodel->referencetemperature)
    7617                                                 +iomodel->latentheat*iomodel->waterfraction[penta_vertex_ids[i]-1];
     7652                                        else nodeinputs[i]=heatcapacity*
     7653                                         (meltingpoint-beta*iomodel->pressure[penta_vertex_ids[i]-1]-referencetemperature)
     7654                                                +latentheat*iomodel->waterfraction[penta_vertex_ids[i]-1];
    76187655                                }
    76197656                                this->inputs->AddInput(new PentaVertexInput(EnthalpyEnum,nodeinputs));

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

    r9320 r9356  
    4545/*}}}*/
    4646/*FUNCTION PentaHook::PentaHook(int in_numanalyses,int matice_id, int matpar_id){{{1*/
    47 PentaHook::PentaHook(int in_numanalyses,int matice_id, int matpar_id){
     47PentaHook::PentaHook(int in_numanalyses,int matice_id, IoModel* iomodel){
     48
     49        /*intermediary: */
     50        int matpar_id;
    4851       
     52        /*retrieve parameters: */
     53        iomodel->parameters->FindParam(&matpar_id,NumberOfElementsEnum); matpar_id++;
     54
    4955        this->numanalyses=in_numanalyses;
    5056        this->hnodes=new Hook*[in_numanalyses];

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

    r4396 r9356  
    88class Hook;
    99class TriaHook;
     10class IoModel;
    1011
    1112class PentaHook{
    … …  
    2021                /*FUNCTION constructors, destructors {{{1*/
    2122                PentaHook();
    22                 PentaHook(int in_numanalyses,int matice_id, int matpar_id);
     23                PentaHook(int in_numanalyses,int matice_id, IoModel* iomodel);
    2324                ~PentaHook();
    2425                void SetHookNodes(int* node_ids,int analysis_counter);

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

    r9320 r9356  
    4141Tria::Tria(int tria_id, int tria_sid, int index, IoModel* iomodel,int nummodels)
    4242        :TriaRef(nummodels)
    43         ,TriaHook(nummodels,index+1,iomodel->numberofelements+1){
     43        ,TriaHook(nummodels,index+1,iomodel){
    4444               
    4545                int i;
    … …  
    534534        /*Retrieve all Inputs and parameters: */
    535535        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    536         this->parameters->FindParam(&artdiff,ArtDiffEnum);
     536        this->parameters->FindParam(&artdiff,ArtificialDiffusivityEnum);
    537537        this->parameters->FindParam(&dim,DimEnum);
    538538        Input* vxaverage_input=NULL;
    … …  
    861861        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    862862        this->parameters->FindParam(&dt,DtEnum);
    863         this->parameters->FindParam(&artdiff,ArtDiffEnum);
     863        this->parameters->FindParam(&artdiff,ArtificialDiffusivityEnum);
    864864        Input* vx_input=inputs->GetInput(HydrologyWaterVxEnum); _assert_(vx_input);
    865865        Input* vy_input=inputs->GetInput(HydrologyWaterVyEnum); _assert_(vy_input);
    … …  
    10201020        this->parameters->FindParam(&dt,DtEnum);
    10211021        this->parameters->FindParam(&dim,DimEnum);
    1022         this->parameters->FindParam(&artdiff,ArtDiffEnum);
     1022        this->parameters->FindParam(&artdiff,ArtificialDiffusivityEnum);
    10231023        Input* vxaverage_input=NULL;
    10241024        Input* vyaverage_input=NULL;
    … …  
    14621462        /*Retrieve all inputs and parameters*/
    14631463        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    1464         this->parameters->FindParam(&num_responses,NumResponsesEnum);
     1464        this->parameters->FindParam(&num_responses,NumCmResponsesEnum);
    14651465        this->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum);
    14661466        Input* thickness_input    = inputs->GetInput(ThicknessEnum);   _assert_(thickness_input);
    … …  
    15281528        /*Retrieve all inputs and parameters*/
    15291529        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    1530         this->parameters->FindParam(&num_responses,NumResponsesEnum);
     1530        this->parameters->FindParam(&num_responses,NumCmResponsesEnum);
    15311531        this->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum);
    15321532        this->parameters->FindParam(&meanvel,MeanvelEnum);
    … …  
    17041704        /*Retrieve all inputs and parameters*/
    17051705        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    1706         this->parameters->FindParam(&num_responses,NumResponsesEnum);
     1706        this->parameters->FindParam(&num_responses,NumCmResponsesEnum);
    17071707        this->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum);
    17081708        this->parameters->FindParam(&meanvel,MeanvelEnum);
    … …  
    27182718        int        *responses = NULL;
    27192719        int         num_responses,resp;
    2720         this->parameters->FindParam(&num_responses,NumResponsesEnum);
     2720        this->parameters->FindParam(&num_responses,NumCmResponsesEnum);
    27212721        this->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum);
    27222722
    … …  
    31123112
    31133113        /*retrieve some parameters: */
    3114         this->parameters->FindParam(&num_controls,NumControlsEnum);
     3114        this->parameters->FindParam(&num_controls,NumControlTypeEnum);
    31153115        this->parameters->FindParam(&control_type,NULL,ControlTypeEnum);
    31163116
    … …  
    32713271        double cmmininputs[3];
    32723272        double cmmaxinputs[3];
     3273        bool   control_analysis=false;
     3274        int    num_control_type;
     3275        double yts;
     3276        int    num_cm_responses;
     3277   
    32733278
    32743279        /*Checks if debuging*/
    … …  
    32773282        /*}}}*/
    32783283
     3284        /*Get parameters: */
     3285        iomodel->parameters->FindParam(&control_analysis,ControlAnalysisEnum);
     3286        iomodel->parameters->FindParam(&num_control_type,NumControlTypeEnum);
     3287        iomodel->parameters->FindParam(&yts,YtsEnum);
     3288        iomodel->parameters->FindParam(&num_cm_responses,NumCmResponsesEnum);
     3289
    32793290        /*Recover vertices ids needed to initialize inputs*/
    32803291        for(i=0;i<3;i++){
    … …  
    32833294
    32843295        /*Control Inputs*/
    3285         if (iomodel->control_analysis && iomodel->control_type){
    3286                 for(i=0;i<iomodel->num_control_type;i++){
     3296        if (control_analysis && iomodel->control_type){
     3297                for(i=0;i<num_control_type;i++){
    32873298                        switch((int)iomodel->control_type[i]){
    32883299                                case DhdtEnum:
    32893300                                        if (iomodel->dhdt){
    3290                                                 for(j=0;j<3;j++)nodeinputs[j]=iomodel->dhdt[tria_vertex_ids[j]-1]/iomodel->yts;
    3291                                                 for(j=0;j<3;j++)cmmininputs[j]=iomodel->cm_min[(tria_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
    3292                                                 for(j=0;j<3;j++)cmmaxinputs[j]=iomodel->cm_max[(tria_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
     3301                                                for(j=0;j<3;j++)nodeinputs[j]=iomodel->dhdt[tria_vertex_ids[j]-1]/yts;
     3302                                                for(j=0;j<3;j++)cmmininputs[j]=iomodel->cm_min[(tria_vertex_ids[j]-1)*num_control_type+i]/yts;
     3303                                                for(j=0;j<3;j++)cmmaxinputs[j]=iomodel->cm_max[(tria_vertex_ids[j]-1)*num_control_type+i]/yts;
    32933304                                                this->inputs->AddInput(new ControlInput(DhdtEnum,TriaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    32943305                                        }
    … …  
    32963307                                case VxEnum:
    32973308                                        if (iomodel->vx){
    3298                                                 for(j=0;j<3;j++)nodeinputs[j]=iomodel->vx[tria_vertex_ids[j]-1]/iomodel->yts;
    3299                                                 for(j=0;j<3;j++)cmmininputs[j]=iomodel->cm_min[(tria_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
    3300                                                 for(j=0;j<3;j++)cmmaxinputs[j]=iomodel->cm_max[(tria_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
     3309                                                for(j=0;j<3;j++)nodeinputs[j]=iomodel->vx[tria_vertex_ids[j]-1]/yts;
     3310                                                for(j=0;j<3;j++)cmmininputs[j]=iomodel->cm_min[(tria_vertex_ids[j]-1)*num_control_type+i]/yts;
     3311                                                for(j=0;j<3;j++)cmmaxinputs[j]=iomodel->cm_max[(tria_vertex_ids[j]-1)*num_control_type+i]/yts;
    33013312                                                this->inputs->AddInput(new ControlInput(VxEnum,TriaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    33023313                                        }
    … …  
    33043315                                case VyEnum:
    33053316                                        if (iomodel->vy){
    3306                                                 for(j=0;j<3;j++)nodeinputs[j]=iomodel->vy[tria_vertex_ids[j]-1]/iomodel->yts;
    3307                                                 for(j=0;j<3;j++)cmmininputs[j]=iomodel->cm_min[(tria_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
    3308                                                 for(j=0;j<3;j++)cmmaxinputs[j]=iomodel->cm_max[(tria_vertex_ids[j]-1)*iomodel->num_control_type+i]/iomodel->yts;
     3317                                                for(j=0;j<3;j++)nodeinputs[j]=iomodel->vy[tria_vertex_ids[j]-1]/yts;
     3318                                                for(j=0;j<3;j++)cmmininputs[j]=iomodel->cm_min[(tria_vertex_ids[j]-1)*num_control_type+i]/yts;
     3319                                                for(j=0;j<3;j++)cmmaxinputs[j]=iomodel->cm_max[(tria_vertex_ids[j]-1)*num_control_type+i]/yts;
    33093320                                                this->inputs->AddInput(new ControlInput(VyEnum,TriaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    33103321                                        }
    … …  
    33133324                                        if (iomodel->drag_coefficient){
    33143325                                                for(j=0;j<3;j++)nodeinputs[j]=iomodel->drag_coefficient[tria_vertex_ids[j]-1];
    3315                                                 for(j=0;j<3;j++)cmmininputs[j]=iomodel->cm_min[(tria_vertex_ids[j]-1)*iomodel->num_control_type+i];
    3316                                                 for(j=0;j<3;j++)cmmaxinputs[j]=iomodel->cm_max[(tria_vertex_ids[j]-1)*iomodel->num_control_type+i];
     3326                                                for(j=0;j<3;j++)cmmininputs[j]=iomodel->cm_min[(tria_vertex_ids[j]-1)*num_control_type+i];
     3327                                                for(j=0;j<3;j++)cmmaxinputs[j]=iomodel->cm_max[(tria_vertex_ids[j]-1)*num_control_type+i];
    33173328                                                this->inputs->AddInput(new ControlInput(DragCoefficientEnum,TriaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    33183329                                        }
    … …  
    33313342                /*Create inputs and add to DataSetInput*/
    33323343                DatasetInput* datasetinput=new DatasetInput(WeightsEnum);
    3333                 for(i=0;i<iomodel->num_cm_responses;i++){
    3334                         for(j=0;j<3;j++)nodeinputs[j]=iomodel->weights[(tria_vertex_ids[j]-1)*iomodel->num_cm_responses+i];
     3344                for(i=0;i<num_cm_responses;i++){
     3345                        for(j=0;j<3;j++)nodeinputs[j]=iomodel->weights[(tria_vertex_ids[j]-1)*num_cm_responses+i];
    33353346                        datasetinput->inputs->AddObject(new TriaVertexInput(WeightsEnum,nodeinputs));
    33363347                }
    … …  
    38933904        double time;
    38943905        TransientInput* transientinput=NULL;
     3906        int    numberofvertices;
     3907        int    numberofelements;
     3908        double yts;
     3909
     3910
     3911        /*Fetch parameters: */
     3912        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     3913        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     3914        iomodel->parameters->FindParam(&yts,YtsEnum);
    38953915
    38963916        /*Branch on type of vector: nodal or elementary: */
    … …  
    39033923
    39043924                /*Are we in transient or static? */
    3905                 if(M==iomodel->numberofvertices){
     3925                if(M==numberofvertices){
    39063926
    39073927                        /*create input values: */
    … …  
    39143934                        this->inputs->AddInput(new TriaVertexInput(vector_enum,nodeinputs));
    39153935                }
    3916                 else if(M==iomodel->numberofvertices+1){
     3936                else if(M==numberofvertices+1){
    39173937                        /*create transient input: */
    39183938                        for(t=0;t<N;t++){ //N is the number of times
    … …  
    39283948
    39293949                                /*time? :*/
    3930                                 time=(double)vector[(M-1)*N+t]*iomodel->yts;
     3950                                time=(double)vector[(M-1)*N+t]*yts;
    39313951
    39323952                                if(t==0) transientinput=new TransientInput(vector_enum);
    … …  
    39353955                        this->inputs->AddInput(transientinput);
    39363956                }
    3937                 else _error_("nodal vector is either numberofnodes (%i), or numberofnodes+1 long. Field provided is %i long",iomodel->numberofvertices,M);
     3957                else _error_("nodal vector is either numberofnodes (%i), or numberofnodes+1 long. Field provided is %i long",numberofvertices,M);
    39383958        }
    39393959        else if(vector_type==2){ //element vector
    39403960                /*Are we in transient or static? */
    3941                 if(M==iomodel->numberofelements){
     3961                if(M==numberofelements){
    39423962
    39433963                        /*static mode: create an input out of the element value: */
    … …  
    52015221        int    tria_type;
    52025222        double nodeinputs[3];
     5223        double yts;
     5224        int    prognostic_DG;
     5225        bool   qmu_analysis;
    52035226
    52045227        /*Checks if debuging*/
    … …  
    52075230        /*}}}*/
    52085231
     5232        /*Fetch parameters: */
     5233        iomodel->parameters->FindParam(&yts,YtsEnum);
     5234        iomodel->parameters->FindParam(&prognostic_DG,PrognosticDGEnum);
     5235        iomodel->parameters->FindParam(&qmu_analysis,QmuAnalysisEnum);
     5236
    52095237        /*Recover element type*/
    5210         if ((analysis_type==PrognosticAnalysisEnum || analysis_type==BalancethicknessAnalysisEnum) && iomodel->prognostic_DG){
     5238        if ((analysis_type==PrognosticAnalysisEnum || analysis_type==BalancethicknessAnalysisEnum) && prognostic_DG){
    52115239                /*P1 Discontinuous Galerkin*/
    52125240                tria_type=P1DGEnum;
    … …  
    52505278                        /*default vx,vy and vz: either observation or 0 */
    52515279                        if(!iomodel->vx){
    5252                                 if (iomodel->vx_obs) for(i=0;i<3;i++)nodeinputs[i]=iomodel->vx_obs[tria_vertex_ids[i]-1]/iomodel->yts;
     5280                                if (iomodel->vx_obs) for(i=0;i<3;i++)nodeinputs[i]=iomodel->vx_obs[tria_vertex_ids[i]-1]/yts;
    52535281                                else                 for(i=0;i<3;i++)nodeinputs[i]=0;
    52545282                                this->inputs->AddInput(new TriaVertexInput(VxEnum,nodeinputs));
    52555283                                this->inputs->AddInput(new TriaVertexInput(VxPicardEnum,nodeinputs));
    5256                                 if(iomodel->qmu_analysis) this->inputs->AddInput(new TriaVertexInput(QmuVxEnum,nodeinputs));
     5284                                if(qmu_analysis) this->inputs->AddInput(new TriaVertexInput(QmuVxEnum,nodeinputs));
    52575285                        }
    52585286                        if(!iomodel->vy){
    5259                                 if (iomodel->vy_obs) for(i=0;i<3;i++)nodeinputs[i]=iomodel->vy_obs[tria_vertex_ids[i]-1]/iomodel->yts;
     5287                                if (iomodel->vy_obs) for(i=0;i<3;i++)nodeinputs[i]=iomodel->vy_obs[tria_vertex_ids[i]-1]/yts;
    52605288                                else                 for(i=0;i<3;i++)nodeinputs[i]=0;
    52615289                                this->inputs->AddInput(new TriaVertexInput(VyEnum,nodeinputs));
    52625290                                this->inputs->AddInput(new TriaVertexInput(VyPicardEnum,nodeinputs));
    5263                                 if(iomodel->qmu_analysis) this->inputs->AddInput(new TriaVertexInput(QmuVyEnum,nodeinputs));
     5291                                if(qmu_analysis) this->inputs->AddInput(new TriaVertexInput(QmuVyEnum,nodeinputs));
    52645292                        }
    52655293                        if(!iomodel->vz){
    5266                                 if (iomodel->vz_obs) for(i=0;i<3;i++)nodeinputs[i]=iomodel->vz_obs[tria_vertex_ids[i]-1]/iomodel->yts;
     5294                                if (iomodel->vz_obs) for(i=0;i<3;i++)nodeinputs[i]=iomodel->vz_obs[tria_vertex_ids[i]-1]/yts;
    52675295                                else                 for(i=0;i<3;i++)nodeinputs[i]=0;
    52685296                                this->inputs->AddInput(new TriaVertexInput(VzEnum,nodeinputs));
    52695297                                this->inputs->AddInput(new TriaVertexInput(VzPicardEnum,nodeinputs));
    5270                                 if(iomodel->qmu_analysis) this->inputs->AddInput(new TriaVertexInput(QmuVzEnum,nodeinputs));
     5298                                if(qmu_analysis) this->inputs->AddInput(new TriaVertexInput(QmuVzEnum,nodeinputs));
    52715299                        }
    52725300                        if(!iomodel->pressure){
    52735301                                for(i=0;i<3;i++)nodeinputs[i]=0;
    5274                                 if(iomodel->qmu_analysis){
     5302                                if(qmu_analysis){
    52755303                                        this->inputs->AddInput(new TriaVertexInput(PressureEnum,nodeinputs));
    52765304                                        this->inputs->AddInput(new TriaVertexInput(QmuPressureEnum,nodeinputs));

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

    r9320 r9356  
    4343/*}}}*/
    4444/*FUNCTION TriaHook::TriaHook(int in_numanalyses,int matice_id, int matpar_id){{{1*/
    45 TriaHook::TriaHook(int in_numanalyses,int matice_id, int matpar_id){
     45TriaHook::TriaHook(int in_numanalyses,int matice_id, IoModel* iomodel){
     46
     47        /*intermediary: */
     48        int matpar_id;
     49
     50        /*retrieve parameters: */
     51        iomodel->parameters->FindParam(&matpar_id,NumberOfElementsEnum); matpar_id++;
    4652       
    4753        this->numanalyses=in_numanalyses;

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

    r4396 r9356  
    77
    88class Hook;
     9class IoModel;
    910
    1011class TriaHook{
    … …  
    1920                /*FUNCTION constructors, destructors {{{1*/
    2021                TriaHook();
    21                 TriaHook(int in_numanalyses,int matice_id, int matpar_id);
     22                TriaHook(int in_numanalyses,int matice_id, IoModel* iomodel);
    2223                ~TriaHook();
    2324                void SetHookNodes(int* node_ids,int analysis_counter);

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

    r9343 r9356  
    1313#include <stdio.h>
    1414#include "./objects.h"
     15#include "./Container/Parameters.h"
    1516#include "../shared/shared.h"
    1617#include "../io/io.h"
    … …  
    1920/*FUNCTION IoModel::IoModel(){{{1*/
    2021IoModel::IoModel(){
     22        this->fid=NULL;
     23        this->parameters=NULL;
     24}
     25/*}}}*/
     26/*FUNCTION IoModel::IoModel(FILE*  iomodel_handle){{{1*/
     27IoModel::IoModel(FILE* iomodel_handle){
     28       
     29        int i,j;
     30
     31        /*First, keep track of the file handle: */
     32        this->fid=iomodel_handle;
     33
     34        /*Initialize and read parameters:*/
     35        this->parameters=new Parameters();
     36        this->FetchParameters(this->parameters); /*this routine goes through the input file, and fetches bools, ints, doubles and strings only, nothing memory intensive*/
     37
     38        /*Initialize the structure: */
    2139        this->IoModelInit();
     40
    2241}
    2342/*}}}*/
    … …  
    4362        xfree((void**)&this->nodeonhutter);
    4463        xfree((void**)&this->nodeonmacayeal);
    45         if (this->dim==3){
    46                 xfree((void**)&this->elements2d);
    47                 xfree((void**)&this->upperelements);
    48                 xfree((void**)&this->lowerelements);
    49                 xfree((void**)&this->nodeonpattyn);
    50         }
     64        xfree((void**)&this->elements2d);
     65        xfree((void**)&this->upperelements);
     66        xfree((void**)&this->lowerelements);
     67        xfree((void**)&this->nodeonpattyn);
    5168        xfree((void**)&this->elementonbed);
    5269        xfree((void**)&this->elementonsurface);
    … …  
    108125
    109126        /*!Delete structure fields: */
    110         xfree((void**)&this->inputfilename);
    111         xfree((void**)&this->outputfilename);
    112         xfree((void**)&this->repository);
    113         xfree((void**)&this->name);
    114127        xfree((void**)&this->riftinfo);
    115128        xfree((void**)&this->penalties);
    … …  
    124137}
    125138/*}}}*/
    126 /*FUNCTION IoModel::IoModel(FILE*  iomodel_handle){{{1*/
    127 IoModel::IoModel(FILE* iomodel_handle){
    128 
    129         /*First, keep track of the fil handle: */
    130         this->fid=iomodel_handle;
    131        
    132         int i,j;
    133                
    134         /*First, initialize the structure: */
    135         this->IoModelInit();
    136        
    137         /*Get all the data that consists of scalars, integers and strings: */
    138 
    139         this->FetchData(&this->name,NameEnum);
    140         this->FetchData(&this->inputfilename,InputfilenameEnum);
    141         this->FetchData(&this->outputfilename,OutputfilenameEnum);
    142         this->FetchData(&this->qmu_analysis,QmuAnalysisEnum);
    143         this->FetchData(&this->control_analysis,ControlAnalysisEnum);
    144         this->FetchData(&this->dim,DimEnum);
    145         /*!Get numberofelements and numberofvertices: */
    146         this->FetchData(&this->numberofvertices,NumberOfNodesEnum);
    147         this->FetchData(&this->numberofelements,NumberOfElementsEnum);
    148         /*!In case we are running 3d, we are going to need the collapsed and non-collapsed 2d meshes, from which the 3d mesh was extruded: */
    149         if (this->dim==3){
    150        
    151                 /*!Deal with 2d mesh: */
    152                 this->FetchData(&this->numberofelements2d,NumberOfElements2DEnum);
    153                 this->FetchData(&this->numberofvertices2d,NumberOfNodes2DEnum);
    154                 this->FetchData(&this->numlayers,NumlayersEnum);
    155         }
    156 
    157 
    158         /*elements type: */
    159         this->FetchData(&this->ishutter,IshutterEnum);
    160         this->FetchData(&this->ismacayealpattyn,IsmacayealpattynEnum);
    161         this->FetchData(&this->isstokes,IsstokesEnum);
    162 
    163         /*!Get drag_type, drag and p,q: */
    164         this->FetchData(&this->drag_type,DragTypeEnum);
    165 
    166         /*!Get materials: */
    167         this->FetchData(&this->rho_water,RhoWaterEnum);
    168         this->FetchData(&this->rho_ice,RhoIceEnum);
    169         this->FetchData(&this->g,GEnum);
    170 
    171         /*Get control parameters: */
    172         this->FetchData(&this->num_control_type,NumControlTypeEnum);
    173         this->FetchData(&this->num_cm_responses,NumCmResponsesEnum);
    174 
    175         /*!Get solution parameters: */
    176         this->FetchData(&this->yts,YtsEnum);
    177         this->FetchData(&this->meanvel,MeanvelEnum);
    178         this->FetchData(&this->epsvel,EpsvelEnum);
    179         this->FetchData(&this->verbose,VerboseBinaryEnum);
    180         this->FetchData(&this->output_frequency,OutputFrequencyEnum);
    181         this->FetchData(&this->prognostic_DG,PrognosticDGEnum);
    182         this->FetchData(&this->nsteps,NstepsEnum);
    183         this->FetchData(&this->eps_cm,EpsCmEnum);
    184         this->FetchData(&this->tolx,TolxEnum);
    185         this->FetchData(&this->cm_gradient,CmGradientEnum);
    186         this->FetchData(&this->eps_res,EpsResEnum);
    187         this->FetchData(&this->eps_rel,EpsRelEnum);
    188         this->FetchData(&this->eps_abs,EpsAbsEnum);
    189         this->FetchData(&this->max_nonlinear_iterations,MaxNonlinearIterationsEnum);
    190         this->FetchData(&this->max_steadystate_iterations,MaxSteadystateIterationsEnum);
    191         this->FetchData(&this->dt,DtEnum);
    192         this->FetchData(&this->ndt,NdtEnum);
    193         this->FetchData(&this->time_adapt,TimeAdaptEnum);
    194         this->FetchData(&this->cfl_coefficient,CflCoefficientEnum);
    195         this->FetchData(&this->hydrostatic_adjustment,HydrostaticAdjustmentEnum);
    196         this->FetchData(&this->penalty_offset,PenaltyOffsetEnum);
    197         this->FetchData(&this->penalty_melting,PenaltyMeltingEnum);
    198         this->FetchData(&this->penalty_lock,PenaltyLockEnum);
    199         this->FetchData(&this->sparsity,SparsityEnum);
    200         this->FetchData(&this->connectivity,ConnectivityEnum);
    201         this->FetchData(&this->lowmem,LowmemEnum);
    202         this->FetchData(&this->viscosity_overshoot,ViscosityOvershootEnum);
    203         this->FetchData(&this->artdiff,ArtificialDiffusivityEnum);
    204         this->FetchData(&this->prognostic_DG,PrognosticDGEnum);
    205         this->FetchData(&this->stokesreconditioning,StokesreconditioningEnum);
    206         this->FetchData(&this->shelf_dampening,ShelfDampeningEnum);
    207         this->FetchData(&this->waitonlock,WaitonlockEnum);
    208         this->FetchData(&this->gl_migration,GlMigrationEnum);
    209         this->FetchData(&this->isdiagnostic,IsdiagnosticEnum);
    210         this->FetchData(&this->isprognostic,IsprognosticEnum);
    211         this->FetchData(&this->isthermal,IsthermalEnum);
    212 
    213         /*!Get thermal parameters: */
    214         this->FetchData(&this->beta,BetaEnum);
    215         this->FetchData(&this->meltingpoint,MeltingpointEnum);
    216         this->FetchData(&this->referencetemperature,ReferencetemperatureEnum);
    217         this->FetchData(&this->latentheat,LatentheatEnum);
    218         this->FetchData(&this->heatcapacity,HeatcapacityEnum);
    219         this->FetchData(&this->thermalconductivity,ThermalconductivityEnum);
    220         this->FetchData(&this->min_thermal_constraints,MinThermalConstraintsEnum);
    221         this->FetchData(&this->min_mechanical_constraints,MinMechanicalConstraintsEnum);
    222         this->FetchData(&this->stabilize_constraints,StabilizeConstraintsEnum);
    223         this->FetchData(&this->mixed_layer_capacity,MixedLayerCapacityEnum);
    224         this->FetchData(&this->thermal_exchange_velocity,ThermalExchangeVelocityEnum);
    225         this->FetchData(&this->basal_melting_rate_correction_apply,BasalMeltingRateCorrectionApplyEnum);
    226         this->FetchData(&this->gl_melting_rate,GlMeltingRateEnum);
    227         this->FetchData(&this->rheology_law,RheologyLawEnum);
    228 
    229         /*!Get hydrology parameters: */         
    230         this->FetchData(&this->hydro_kn,HydroKnEnum);                   
    231         this->FetchData(&this->hydro_p,HydroPEnum);             
    232         this->FetchData(&this->hydro_q,HydroQEnum);             
    233         this->FetchData(&this->hydro_CR,HydroCREnum);           
    234         this->FetchData(&this->hydro_n,HydroNEnum);
    235 
    236         /*qmu: */
    237         if(this->qmu_analysis){
    238                 this->FetchData(&this->numberofvariables,NumberOfVariablesEnum);
    239                 this->FetchData(&this->numberofresponses,NumberOfResponsesEnum);
    240                 this->FetchData(&this->qmu_npart,NpartEnum);
    241                 this->FetchData(&this->qmu_save_femmodel,QmuSaveFemmodelEnum);
    242         }
    243 
    244         /*i/o: */
    245         this->FetchData(&this->io_gather,IoGatherEnum);
    246 }
    247 /*}}}*/
    248139/*FUNCTION IoModel::IoModelInit{{{1*/
    249140void IoModel::IoModelInit(void){
     141
     142        /*exterior data: */
     143        this->nodecounter=0;
     144        this->loadcounter=0;
     145        this->constraintcounter=0;
    250146       
    251147        /*!initialize all pointers to 0: */
    252         this->name=NULL;
    253         this->inputfilename=NULL;
    254         this->outputfilename=NULL;
    255         this->repository=NULL;
    256         this->qmu_analysis=0;
    257         this->control_analysis=0;
    258         this->numberofvariables=0;
    259         this->numvariabledescriptors=0;
    260         this->numberofresponses=0;
    261         this->numresponsedescriptors=0;
    262         this->qmu_npart=0;
    263         this->numberofelements=0;
    264         this->numberofvertices=0;
    265148        this->x=NULL;
    266149        this->y=NULL;
    … …  
    270153        this->elements_type=NULL;
    271154        this->vertices_type=NULL;
    272         this->numberofvertices2d=0;
    273155        this->elements2d=NULL;
    274         this->numlayers=0;
    275156        this->upperelements=NULL;
    276157        this->lowerelements=NULL;
    … …  
    278159        this->nodeonmacayeal=NULL;
    279160        this->nodeonpattyn=NULL;
    280         this->io_gather=1;
    281161       
    282162        this->vx_obs=NULL;
    … …  
    289169        this->temperature=NULL;
    290170        this->waterfraction=NULL;
    291         this->gl_melting_rate=0;
    292171        this->basal_melting_rate=NULL;
    293172        this->watercolumn=NULL;
    294173        this->basal_melting_rate_correction=NULL;
    295         this->basal_melting_rate_correction_apply=0;
    296174        this->geothermalflux=NULL;
    297175        this->elementonbed=NULL;
    … …  
    313191        this->nodeonwater=NULL;
    314192
    315         this->drag_type=0;
    316193        this->drag_coefficient=NULL;
    317194        this->drag_p=NULL;
    … …  
    319196       
    320197       
    321         this->numberofpressureloads=0;
    322198        this->pressureload=NULL;
    323199        this-> spcvx=NULL;
    … …  
    328204        this-> spcwatercolumn=NULL;
    329205        this-> diagnostic_ref=NULL;
    330         this->numberofedges=0;
    331206        this->edges=NULL;
    332207       
    333208        /*!materials: */
    334         this->rho_water=0;
    335         this->rho_ice=0;
    336         this->g=0;
    337209        this->rheology_n=NULL;
    338210        this->rheology_B=NULL;
    339         this->rheology_law=0;
    340211
    341212        /*!solution parameters: */
    … …  
    344215        this->weights=NULL;
    345216        this->cm_jump=NULL;
    346         this->meanvel=0;
    347         this->epsvel=0;
    348         this->nsteps=0;
    349         this->eps_cm=0;
    350         this->tolx=0;
    351217        this->maxiter=NULL;
    352218        this->cm_min=NULL;
    353219        this->cm_max=NULL;
    354         this->cm_gradient=0;
    355         this->verbose=0;
    356         this->output_frequency=0;
    357         this->eps_res=0;
    358         this->eps_rel=0;
    359         this->eps_abs=0;
    360         this->max_nonlinear_iterations=0;
    361         this->max_steadystate_iterations=0;
    362         this->dt=0;
    363         this->ndt=0;
    364         this->time_adapt=0;
    365         this->cfl_coefficient=0;
    366         this->hydrostatic_adjustment=0;
    367         this->gl_migration=0;
    368         this->penalty_offset=0;
    369         this->penalty_melting=0;
    370         this->penalty_lock=0;
    371         this->sparsity=0;
    372         this->connectivity=0;
    373         this->lowmem=0;
    374220        this->optscal=NULL;
    375         this->yts=0;
    376         this->viscosity_overshoot=0;
    377         this->artdiff=0;
    378         this->prognostic_DG=0;
    379         this->stokesreconditioning=0;
    380         this->shelf_dampening=0;
    381         this->waitonlock=0;
    382         this->isdiagnostic=0;
    383         this->isprognostic=0;
    384         this->isthermal=0;
    385 
    386         /*!thermal parameters: */
    387         this->beta=0;
    388         this->meltingpoint=0;
    389         this->referencetemperature=0;
    390         this->latentheat=0;
    391         this->heatcapacity=0;
    392         this->thermalconductivity=0;
    393         this->min_thermal_constraints=0;
    394         this->min_mechanical_constraints=0;
    395         this->stabilize_constraints=0;
    396         this->mixed_layer_capacity=0;
    397         this->thermal_exchange_velocity=0;
    398 
    399        
    400         this->numrifts=0;
     221
    401222        this->riftinfo=NULL;
    402223
    403224        /*!penalties: */
    404         this->numpenalties=0;
    405225        this->penalties=NULL;
    406226
    … …  
    410230        this->surface_ablation_rate=NULL;
    411231        this->dhdt=NULL;
    412 
    413         /*elements type: */
    414         this->ishutter=0;
    415         this->ismacayealpattyn=0;
    416         this->isstokes=0;
    417232
    418233        /*exterior data: */
    … …  
    422237        this->singlenodetoelementconnectivity=NULL;
    423238        this->numbernodetoelementconnectivity=NULL;
    424         this->nodecounter=0;
    425         this->loadcounter=0;
    426         this->constraintcounter=0;
    427239}
    428240/*}}}*/
    … …  
    435247
    436248        /*output: */
    437         bool  boolean;
     249        int   booleanint;
    438250        int   code;
    439251       
    440252        /*Set file pointer to beginning of the data: */
    441253        fid=this->SetFilePointerToData(&code,NULL,data_enum);
    442        
     254
    443255        if(code!=1)_error_("%s%s","IoModel::FetchData expecting a boolean for enum ",EnumToStringx(data_enum));
    444256       
    445257        /*We have to read a boolean from disk. */
    446258        if(my_rank==0){ 
    447                 if(fread(&boolean,sizeof(bool),1,fid)!=1) _error_(" could not read boolean ");
    448         }
    449 
    450         MPI_Bcast(&boolean,1,MPI_BYTE,0,MPI_COMM_WORLD);
    451 
     259                if(fread(&booleanint,sizeof(int),1,fid)!=1) _error_(" could not read boolean ");
     260        }
     261        MPI_Bcast(&booleanint,1,MPI_INT,0,MPI_COMM_WORLD);
     262
     263        /*cast to bool: */
    452264        /*Assign output pointers: */
    453         *pboolean=boolean;
     265        *pboolean=(bool)booleanint;
    454266
    455267}
    … …  
    558370}
    559371/*}}}*/
    560 /*FUNCTION IoModel::FetchData(double**  pdoublematrix,int* pM,int* pN,int data_enum){{{1*/
    561 void  IoModel::FetchData(double** pmatrix,int* pM,int* pN,int data_enum){
     372/*FUNCTION IoModel::FetchData(int**  pintegerematrix,int* pM,int* pN,int data_enum){{{1*/
     373void  IoModel::FetchData(int** pmatrix,int* pM,int* pN,int data_enum){
    562374
    563375        extern int my_rank;
    564376        extern int num_procs;
     377        int i,j;
    565378
    566379        /*output: */
    567380        int M,N;
    568381        double* matrix=NULL;
     382        int*    integer_matrix=NULL;
    569383        int code=0;
    570384        int vector_type=0;
    571        
    572385       
    573386       
    … …  
    604417        }
    605418
     419        /*Now cast to integer: */
     420        if(M*N){
     421                integer_matrix=(int*)xmalloc(M*N*sizeof(int));
     422                for (i=0;i<M;i++){
     423                        for (j=0;j<N;j++){
     424                                integer_matrix[i*N+j]=(int)matrix[i*N+j];
     425                        }
     426                }
     427        }
     428        else{
     429                integer_matrix=NULL;
     430        }
     431        /*Free ressources:*/
     432        xfree((void**)&matrix);
     433
     434        /*Assign output pointers: */
     435        *pmatrix=integer_matrix;
     436        if (pM)*pM=M;
     437        if (pN)*pN=N;
     438
     439}
     440/*}}}*/
     441/*FUNCTION IoModel::FetchData(double**  pdoublematrix,int* pM,int* pN,int data_enum){{{1*/
     442void  IoModel::FetchData(double** pmatrix,int* pM,int* pN,int data_enum){
     443
     444        extern int my_rank;
     445        extern int num_procs;
     446
     447        /*output: */
     448        int M,N;
     449        double* matrix=NULL;
     450        int code=0;
     451        int vector_type=0;
     452       
     453       
     454       
     455        /*Set file pointer to beginning of the data: */
     456        fid=this->SetFilePointerToData(&code,&vector_type,data_enum);
     457
     458        if((code!=5) && (code!=6) && (code!=7))_error_("%s%s","IoModel::FetchData expecting a double, integer or boolean matrix for enum ",EnumToStringx(data_enum));
     459       
     460        /*Now fetch: */
     461
     462        /*We have to read a matrix from disk. First read the dimensions of the matrix, then the whole matrix: */
     463        /*numberofelements: */
     464        if(my_rank==0){ 
     465                if(fread(&M,sizeof(int),1,fid)!=1) _error_("could not read number of rows for matrix ");
     466        }
     467
     468        MPI_Bcast(&M,1,MPI_INT,0,MPI_COMM_WORLD);
     469
     470        if(my_rank==0){ 
     471                if(fread(&N,sizeof(int),1,fid)!=1) _error_("could not read number of columns for matrix ");
     472        }
     473        MPI_Bcast(&N,1,MPI_INT,0,MPI_COMM_WORLD);
     474
     475        /*Now allocate matrix: */
     476        if(M*N){
     477                matrix=(double*)xmalloc(M*N*sizeof(double));
     478
     479                /*Read matrix on node 0, then broadcast: */
     480                if(my_rank==0){ 
     481                        if(fread(matrix,M*N*sizeof(double),1,fid)!=1) _error_("could not read matrix ");
     482                }
     483               
     484                MPI_Bcast(matrix,M*N,MPI_DOUBLE,0,MPI_COMM_WORLD);
     485        }
     486
    606487        /*Assign output pointers: */
    607488        *pmatrix=matrix;
    … …  
    757638        *pndims=ndims;
    758639        *pnumrecords=numrecords;
     640}
     641/*}}}*/
     642/*FUNCTION IoModel::FetchParameters(Parameters* parameters){{{1*/
     643void  IoModel::FetchParameters(Parameters* parameters){
     644
     645        extern int my_rank;
     646        extern int num_procs;
     647       
     648        /*record descriptions; */
     649        int record_enum;
     650        int record_length;
     651        int record_code; //1 to 7 number
     652
     653        /*records: */
     654        int  booleanint=0;
     655        int  integer=0;
     656        double scalar=0;
     657        char* string=NULL;
     658        int   string_size;
     659
     660 
     661        /*Go find in the binary file, the position of the data we want to fetch: */
     662        if(my_rank==0){ //cpu 0{{{2
     663       
     664                /*First set FILE* position to the beginning of the file: */
     665                fseek(this->fid,0,SEEK_SET);
     666
     667                /*Now march through file looking for the correct data identifiers (bool,int,double or string): */
     668                for(;;){
     669                        if(fread(&record_enum,sizeof(int),1,this->fid)==0){
     670
     671                                /*Ok, we have reached the end of the file. break: */
     672                                record_code=0; //0 means bailout
     673                                MPI_Bcast(&record_code,1,MPI_INT,0,MPI_COMM_WORLD);  /*tell others cpus we are bailing: */
     674                                break;
     675                        }
     676                        else{
     677                       
     678                                /* Read the record length and the data type code: */
     679                                fread(&record_length,sizeof(int),1,this->fid);
     680                                fread(&record_code,sizeof(int),1,this->fid);
     681                                       
     682                                /*Tell other cpus what we are doing: */
     683                                MPI_Bcast(&record_code,1,MPI_INT,0,MPI_COMM_WORLD);  /*tell other cpus what we are going to do: */
     684
     685                                /*Tell other cpus the name of the data, then branch according to the data type: */
     686                                MPI_Bcast(&record_enum,1,MPI_INT,0,MPI_COMM_WORLD); 
     687                                MPI_Bcast(&record_length,1,MPI_INT,0,MPI_COMM_WORLD); 
     688                               
     689
     690                                switch(record_code){
     691                                        case 1:
     692                                                /*Read the boolean and broadcast it to other cpus:*/
     693                                                if(fread(&booleanint,sizeof(int),1,this->fid)!=1) _error_(" could not read boolean ");
     694                                                MPI_Bcast(&booleanint,1,MPI_INT,0,MPI_COMM_WORLD);
     695
     696                                                /*create BoolParam: */
     697                                                parameters->AddObject(new BoolParam(record_enum,(bool)booleanint)); //cast to boolean
     698
     699                                                break;
     700                                        case 2:
     701                                                /*Read the integer and broadcast it to other cpus:*/
     702                                                if(fread(&integer,sizeof(int),1,this->fid)!=1) _error_(" could not read integer ");
     703                                                MPI_Bcast(&integer,1,MPI_INT,0,MPI_COMM_WORLD);
     704
     705                                                /*create IntParam: */
     706                                                parameters->AddObject(new IntParam(record_enum,integer));
     707
     708                                                break;
     709                                        case 3:
     710                                                /*Read the scalar and broadcast it to other cpus:*/
     711                                                if(fread(&scalar,sizeof(double),1,this->fid)!=1) _error_(" could not read scalar ");
     712                                                MPI_Bcast(&scalar,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
     713
     714                                                /*create DoubleParam: */
     715                                                parameters->AddObject(new DoubleParam(record_enum,scalar));
     716
     717                                                break;
     718                                        case 4:
     719                                                /*We have to read a string from disk. First read the dimensions of the string, then the string: */
     720                                                if(fread(&string_size,sizeof(int),1,this->fid)!=1) _error_(" could not read length of string ");
     721                                                MPI_Bcast(&string_size,1,MPI_INT,0,MPI_COMM_WORLD);
     722
     723                                                if(string_size){
     724                                                        string=(char*)xmalloc((string_size+1)*sizeof(char));
     725                                                        string[string_size]='\0';
     726
     727                                                        /*Read string, then broadcast: */
     728                                                        if(fread(string,string_size*sizeof(char),1,this->fid)!=1)_error_("  could not read string ");
     729                                                        MPI_Bcast(string,string_size,MPI_CHAR,0,MPI_COMM_WORLD);
     730                                                }
     731                                                else{
     732                                                        string=(char*)xmalloc(sizeof(char));
     733                                                        string[0]='\0';
     734                                                }
     735                                               
     736                                                /*Add string to parameters: */
     737                                                parameters->AddObject(new StringParam(record_enum,string));
     738
     739                                                break;
     740                                        case 5:
     741                                                        /*We are not interested in this record, too memory intensive. Skip it: */
     742                                                        /*skip: */
     743                                                        fseek(fid,-sizeof(int),SEEK_CUR); //backtrak 1 integer
     744                                                        fseek(fid,record_length,SEEK_CUR);
     745                                                        break;
     746                                        case 6:
     747                                                        /*We are not interested in this record, too memory intensive. Skip it: */
     748                                                        /*skip: */
     749                                                        fseek(fid,-sizeof(int),SEEK_CUR); //backtrak 1 integer
     750                                                        fseek(fid,record_length,SEEK_CUR);
     751                                                        break;
     752                                        case 7:
     753                                                        /*We are not interested in this record, too memory intensive. Skip it: */
     754                                                        /*skip: */
     755                                                        fseek(fid,-sizeof(int),SEEK_CUR); //backtrak 1 integer
     756                                                        fseek(fid,record_length,SEEK_CUR);
     757                                                        break;
     758
     759                                        case 8:
     760                                                        /*We are not interested in this record, too memory intensive. Skip it: */
     761                                                        /*skip: */
     762                                                        fseek(fid,-sizeof(int),SEEK_CUR); //backtrak 1 integer
     763                                                        fseek(fid,record_length,SEEK_CUR);
     764                                                        break;
     765
     766                                        case 9:
     767                                                        /*We are not interested in this record, too memory intensive. Skip it: */
     768                                                        /*skip: */
     769                                                        fseek(fid,-sizeof(int),SEEK_CUR); //backtrak 1 integer
     770                                                        fseek(fid,record_length,SEEK_CUR);
     771                                                        break;
     772
     773                                        default:
     774                                                _error_("%s%i","unknown record type:",record_code);
     775                                                break;;
     776                                }
     777                        }
     778                }
     779        } //}}}
     780        else{ //cpu ~0 {{{2
     781                for(;;){ //wait on cpu 0
     782                        MPI_Bcast(&record_code,1,MPI_INT,0,MPI_COMM_WORLD);  /*get from cpu 0 what we are going to do: */
     783                        if(record_code==0){
     784                                break; //we are done, break from the loop
     785                        }
     786                        else{
     787                                MPI_Bcast(&record_enum,1,MPI_INT,0,MPI_COMM_WORLD);   //get from cpu 0 name of the data
     788                                MPI_Bcast(&record_length,1,MPI_INT,0,MPI_COMM_WORLD); 
     789                                switch(record_code){
     790                                case 1:
     791                                        /*boolean. get it from cpu 0 */
     792                                        MPI_Bcast(&booleanint,1,MPI_INT,0,MPI_COMM_WORLD);
     793                                               
     794                                        /*create BoolParam: */
     795                                        parameters->AddObject(new BoolParam(record_enum,(bool)booleanint)); //cast to a boolean
     796                                        break;
     797
     798                                case 2:
     799                                        /*integer. get it from cpu 0 */
     800                                        MPI_Bcast(&integer,1,MPI_INT,0,MPI_COMM_WORLD);
     801                                               
     802                                        /*create IntParam: */
     803                                        parameters->AddObject(new IntParam(record_enum,integer));
     804
     805                                        break;
     806                                case 3:
     807                                        /*scalar. get it from cpu 0 */
     808                                        MPI_Bcast(&scalar,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
     809                                               
     810                                        /*create DoubleParam: */
     811                                        parameters->AddObject(new DoubleParam(record_enum,scalar));
     812
     813                                        break;
     814                                case 4:
     815                                        MPI_Bcast(&string_size,1,MPI_INT,0,MPI_COMM_WORLD);
     816                                        if(string_size){
     817                                                string=(char*)xmalloc((string_size+1)*sizeof(char));
     818                                                string[string_size]='\0';
     819
     820                                                /*Read string from cpu 0: */
     821                                                MPI_Bcast(string,string_size,MPI_CHAR,0,MPI_COMM_WORLD);
     822                                        }
     823                                        else{
     824                                                string=(char*)xmalloc(sizeof(char));
     825                                                string[0]='\0';
     826                                        }
     827                                        /*Add string to parameters: */
     828                                        parameters->AddObject(new StringParam(record_enum,string));
     829
     830                                        break;
     831                                case 5: break; //do nothing. not interested in this type of data, which is memory intensive.
     832                                case 6: break; //do nothing. not interested in this type of data, which is memory intensive.
     833                                case 7: break; //do nothing. not interested in this type of data, which is memory intensive.
     834                                case 8: break; //do nothing. not interested in this type of data, which is memory intensive.
     835                                case 9: break; //do nothing. not interested in this type of data, which is memory intensive.
     836
     837                                default:
     838                                        _error_("%s%i","unknown record type:",record_code);
     839                                        break;;
     840                                }
     841
     842
     843                        }
     844                }
     845        } //}}}
    759846}
    760847/*}}}*/
    … …  
    821908}
    822909/*}}}*/
    823 /*FUNCTION IoModel::Echo{{{1*/
    824 void IoModel::Echo(int which_part,int rank) {
    825 
    826         //which_part  determines what gets echoed, otherwise, we'll get too much output.
    827         //1-> penalties
    828 
    829         int i,j;
    830 
    831         if(which_part==1 && my_rank==rank && this->dim==3){
    832                 printf("IoModel penalties: \n");
    833                 printf("   number of penalties: %i\n",this->numpenalties);
    834                 printf("   nodes: \n");
    835 
    836                 for(i=0;i<this->numpenalties;i++){
    837                         for(j=0;j<this->numlayers;j++){
    838                                 printf("%i ",(int)*(this->penalties+this->numlayers*i+j));
    839                         }
    840                         printf("\n");
    841                 }
    842         }
    843 }
    844 /*}}}*/
    845910/*FUNCTION IoModel::FetchDataToInput(Elements* elements,IoModel* iomodel,int vector_enum,int default_vector_enum,double default_value{{{1*/
    846911void IoModel::FetchDataToInput(Elements* elements,int vector_enum,int default_vector_enum,double default_value){
    … …  
    858923        int     nods;
    859924        int     nel;
     925        int     numberofelements;
    860926
    861927
    … …  
    868934        int     M,N;
    869935
     936        /*Fetch parameters: */
     937        this->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     938
    870939        /*First of, find the record for the enum, and get code  of data type: */
    871940        fid=this->SetFilePointerToData(&code, &vector_layout,vector_enum);
    … …  
    877946                        /*Add boolean constant input to all elements: */
    878947                        counter=0;
    879                         for (i=0;i<this->numberofelements;i++){
     948                        for (i=0;i<numberofelements;i++){
    880949                                if(this->my_elements[i]){
    881950                                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    890959                        /*Add integer constant input to all elements: */
    891960                        counter=0;
    892                         for (i=0;i<this->numberofelements;i++){
     961                        for (i=0;i<numberofelements;i++){
    893962                                if(this->my_elements[i]){
    894963                                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    903972                        /*Add double constant input to all elements: */
    904973                        counter=0;
    905                         for (i=0;i<this->numberofelements;i++){
     974                        for (i=0;i<numberofelements;i++){
    906975                                if(this->my_elements[i]){
    907976                                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    9401009                        /*Create inputs:*/
    9411010                        counter=0;
    942                         for (i=0;i<this->numberofelements;i++){
     1011                        for (i=0;i<numberofelements;i++){
    9431012                                if(this->my_elements[i]){
    9441013                                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    9781047                        /*Create inputs:*/
    9791048                        counter=0;
    980                         for (i=0;i<this->numberofelements;i++){
     1049                        for (i=0;i<numberofelements;i++){
    9811050                                if(this->my_elements[i]){
    9821051                                        Element* element=(Element*)elements->GetObjectByOffset(counter);
    … …  
    10161085                        /*Create inputs:*/
    10171086                        counter=0;
    1018                         for (i=0;i<this->numberofelements;i++){
     1087                        for (i=0;i<numberofelements;i++){
    10191088                                if(this->my_elements[i]){
    10201089                                        Element* element=(Element*)elements->GetObjectByOffset(counter);

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

    r9343 r9356  
    1616        private:
    1717                FILE* fid; //pointer to input file
     18
    1819        public:
     20                Parameters* parameters; //this dataset holds all double, int, bool and char* parameters read in from the input file.*
    1921                /*Data: {{{1*/
    20                 char*   name;
    21                 char*   inputfilename;
    22                 char*   outputfilename;
    23                 char*   repository;
    24                 int     dim;
    25                 int     qmu_analysis;
    26                 int     control_analysis;
     22
     23                int     nodecounter; //keep track of how many nodes are being created in each analysis type
     24                int     loadcounter; //keep track of how many loads are being created in each analysis type
     25                int     constraintcounter; //keep track of how many constraints are being created in each analysis type
    2726
    2827                /*2d mesh: */
    29                 int     numberofelements;
    30                 int     numberofvertices;
    3128                double* x;
    3229                double* y;
    … …  
    3835
    3936                /*3d mesh: */
    40                 int     numberofvertices2d;
    41                 int     numberofelements2d;
    4237                double* elements2d;
    43                 int     numlayers;
    4438                double* upperelements;
    4539                double* lowerelements;
    4640
    4741                /*elements type: */
    48                 int     ishutter;
    49                 int     ismacayealpattyn;
    50                 int     isstokes;
    5142                double* nodeonhutter;
    5243                double* nodeonmacayeal;
    … …  
    6152                double* waterfraction;
    6253
    63                 /*i/o: */
    64                 int    io_gather;
    65 
    6654                /*observations: */
    6755                double*  thickness_obs;
    … …  
    6957                double*  vy_obs;
    7058                double*  vz_obs;
    71 
    72                 /*qmu: */
    73                 int      numberofvariables;
    74                 int      numvariabledescriptors;
    75                 int      numberofresponses;
    76                 int      numresponsedescriptors;
    77                 int      qmu_npart;
    7859
    7960                /*geometry: */
    … …  
    9677
    9778                /*friction: */
    98                 int     drag_type;
    9979                double* drag_coefficient;
    10080                double* drag_p;
    … …  
    10282
    10383                /*boundary conditions: */
    104                 int     numberofpressureloads;
    10584                double* pressureload;
    10685                double* spcvx;
    … …  
    11291                double* diagnostic_ref;
    11392                double* geothermalflux;
    114                 int     numberofedges;
    11593                double* edges;
    11694
    11795                /*materials: */
    118                 double  rho_water,rho_ice;
    119                 double  g;
    12096                double* rheology_B;
    12197                double* rheology_n;
    122                 int     rheology_law;
    12398
    12499                /*numerical parameters: */
    125                 double  meanvel;
    126                 double  epsvel;
    127                 int     artdiff;
    128                 int     prognostic_DG;
    129                 double  viscosity_overshoot;
    130                 double  stokesreconditioning;
    131                 int     shelf_dampening;
    132100                double* cm_min;
    133101                double* cm_max;
    134                 int     cm_gradient;;
    135102
    136103                /*control methods: */
    137                 int      num_cm_responses;
    138                 int      num_control_type;
    139104                double*  control_type;
    140105
    … …  
    143108                double* weights;
    144109                double* cm_jump;
    145                 int     nsteps;
    146                 double  eps_cm;
    147                 double  tolx;
    148110                double* maxiter;
    149                 int     verbose;
    150                 int     output_frequency;
    151                 double  eps_res;
    152                 double  eps_rel;
    153                 double  eps_abs;
    154                 double  max_nonlinear_iterations;
    155                 double  max_steadystate_iterations;
    156                 double  dt,ndt;
    157                 int     time_adapt;
    158                 double  cfl_coefficient;
    159                 int     hydrostatic_adjustment;
    160                 double  penalty_offset;
    161                 double  penalty_melting;
    162                 int     penalty_lock;
    163                 double  sparsity;
    164                 int     connectivity;
    165                 int     lowmem;
    166111                double* optscal;
    167                 double  yts;
    168                 double  waitonlock;
    169                 int     isdiagnostic;
    170                 int     isprognostic;
    171                 int     isthermal;
    172 
    173                 /*thermal parameters: */
    174                 double beta;
    175                 double meltingpoint;
    176                 double referencetemperature;
    177                 double latentheat;
    178                 double  heatcapacity,thermalconductivity;
    179                 int    min_thermal_constraints;
    180                 int    min_mechanical_constraints;
    181                 int    stabilize_constraints;
    182                 double mixed_layer_capacity;
    183                 double thermal_exchange_velocity;
    184112
    185113                /*rifts: */
    186                 int      numrifts;
    187114                double*  riftinfo;
    188115
    189116                /*penalties: */
    190                 int      numpenalties;
    191117                double*  penalties;
    192118
    … …  
    194120                double*  basal_melting_rate;
    195121                double*  watercolumn;
    196                 double   gl_melting_rate;
    197122                double*  basal_melting_rate_correction;
    198                 int      basal_melting_rate_correction_apply;
    199123                double*  surface_mass_balance;
    200124                double*  surface_accumulation_rate;
    … …  
    202126
    203127                double*  dhdt;
    204 
    205                 /*qmu: */
    206                 int      qmu_save_femmodel;
    207 
    208                 /*grounding line migration: */
    209                 int      gl_migration;
    210128
    211129                /*exterior partitioning data, to be carried around: */
    … …  
    215133                int*    singlenodetoelementconnectivity;
    216134                int*    numbernodetoelementconnectivity;
    217                 int     nodecounter; //keep track of how many nodes are being created in each analysis type
    218                 int     loadcounter; //keep track of how many loads are being created in each analysis type
    219                 int     constraintcounter; //keep track of how many constraints are being created in each analysis type
    220135
    221                 /*hydrology parameters: */               
    222                 double hydro_kn;                 
    223                 double hydro_p;         
    224                 double hydro_q;         
    225                 double hydro_CR;                         
    226                 double hydro_n;
    227136                 /*}}}*/
    228137                /*Methods: {{{1*/
    … …  
    231140                IoModel(FILE* iomodel_handle);
    232141                void IoModelInit(void);
    233                 void Echo(int which_part,int rank);
    234142                /*}}}*/
    235143                /*Input/Output:{{{1*/
    … …  
    238146                void  FetchData(double*   pscalar,int data_enum);
    239147                void  FetchData(char**    pstring,int data_enum);
     148                void  FetchData(int** pmatrix,int* pM,int* pN,int data_enum);
    240149                void  FetchData(double**  pscalarmatrix,int* pM,int* pN,int data_enum);
    241150                void  FetchData(char***   pstringarray,int* pnumstrings,int data_enum);
    242151                void  FetchData(double*** pmatrixarray,int** pmdims,int** pndims, int* pnumrecords,int data_enum);
     152                void  FetchParameters(Parameters* parameters);
    243153                FILE* SetFilePointerToData(int* pcode,int* pvector_type, int data_enum);
    244154                void  FetchDataToInput(Elements* elements,int vector_enum,int default_vector_enum=NoneEnum,double default_value=0);

  • issm/trunk/src/c/objects/Loads/Icefront.cpp

    r9320 r9356  
    4343        int segment_width;
    4444        int element;
    45         int num_nodes;
     45        int num_nodes;
     46        int dim;
     47        int numberofelements;
    4648
    4749        /*icefront constructor data: */
    … …  
    5052        int  icefront_node_ids[NUMVERTICESQUA]; //initialize with largest size
    5153        int  icefront_fill;
     54       
     55        /*find parameters: */
     56        iomodel->parameters->FindParam(&dim,DimEnum);
     57        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
    5258
    5359        /*First, retrieve element index and element type: */
    54         if (iomodel->dim==2){
     60        if (dim==2){
    5561                segment_width=4;
    5662        }
    … …  
    6268        /*Build ids for hook constructors: */
    6369        icefront_eid=(int) *(iomodel->pressureload+segment_width*i+segment_width-2); //matlab indexing
    64         icefront_mparid=iomodel->numberofelements+1; //matlab indexing
     70        icefront_mparid=numberofelements+1; //matlab indexing
    6571
    6672        if (in_icefront_type==MacAyeal2dIceFrontEnum || in_icefront_type==MacAyeal3dIceFrontEnum){

  • issm/trunk/src/c/objects/Loads/Numericalflux.cpp

    r9320 r9356  
    5151        int   numericalflux_type;
    5252
     53        int    numberofelements;
     54
     55        /*Fetch parameters: */
     56        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     57
    5358        /* Get MatPar id */
    54         numericalflux_mparid=iomodel->numberofelements+1; //matlab indexing
     59        numericalflux_mparid=numberofelements+1; //matlab indexing
    5560
    5661        /*First, see wether this is an internal or boundary edge (if e2=NaN)*/

  • issm/trunk/src/c/objects/Loads/Pengrid.cpp

    r9320 r9356  
    4949        int pengrid_element_id;
    5050
     51        int numberofvertices;
     52        int numberofelements;
     53
     54        /*Fetch parameters: */
     55        iomodel->parameters->FindParam(&numberofvertices,NumberOfVerticesEnum);
     56        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     57
    5158        /*Some checks if debugging activated*/
    5259        _assert_(iomodel->singlenodetoelementconnectivity);
    53         _assert_(index>=0 && index<iomodel->numberofvertices);
     60        _assert_(index>=0 && index<numberofvertices);
    5461        _assert_(id);
    5562
    … …  
    6269        pengrid_element_id=iomodel->singlenodetoelementconnectivity[index];
    6370        _assert_(pengrid_element_id);
    64         pengrid_matpar_id=iomodel->numberofelements+1; //refers to the constant material parameters object
     71        pengrid_matpar_id=numberofelements+1; //refers to the constant material parameters object
    6572
    6673        this->hnode=new Hook(&pengrid_node_id,1);

  • issm/trunk/src/c/objects/Loads/Riftfront.cpp

    r9320 r9356  
    4848        double riftfront_fractionincrement;
    4949        bool   riftfront_shelf;
     50        int    numberofelements;
     51        int    penalty_lock;
    5052
    5153        /*intermediary: */
    5254        int el1    ,el2;
    5355        int node1  ,node2;
     56
     57        /*Fetch parameters: */
     58        iomodel->parameters->FindParam(&numberofelements,NumberOfElementsEnum);
     59        iomodel->parameters->FindParam(&penalty_lock,PenaltyLockEnum);
    5460
    5561        /*Ok, retrieve all the data needed to add a penalty between the two nodes: */
    … …  
    6975        riftfront_elem_ids[0]=el1;
    7076        riftfront_elem_ids[1]=el2;
    71         riftfront_matpar_id=iomodel->numberofelements+1; //matlab indexing
     77        riftfront_matpar_id=numberofelements+1; //matlab indexing
    7278
    7379        /*Hooks: */
    … …  
    8187        this->counter=0;
    8288        this->prestable=0;
    83         this->penalty_lock=iomodel->penalty_lock;
     89        this->penalty_lock=penalty_lock;
    8490        this->material_converged=0;
    8591        this->normal[0]=*(iomodel->riftinfo+RIFTINFOSIZE*i+4);

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

    r9320 r9356  
    655655        int i,j;
    656656
     657        int    dim;
     658        int    control_analysis;
     659        int    num_control_type;
     660
     661        /*Fetch parameters: */
     662        iomodel->parameters->FindParam(&dim,DimEnum);
     663        iomodel->parameters->FindParam(&control_analysis,ControlAnalysisEnum);
     664        iomodel->parameters->FindParam(&num_control_type,NumControlTypeEnum);
     665
     666
    657667        /*if 2d*/
    658         if(iomodel->dim==2){
     668        if(dim==2){
    659669
    660670                /*Intermediaries*/
    … …  
    677687
    678688                /*Control Inputs*/
    679                 if (iomodel->control_analysis && iomodel->control_type){
    680                         for(i=0;i<iomodel->num_control_type;i++){
     689                if (control_analysis && iomodel->control_type){
     690                        for(i=0;i<num_control_type;i++){
    681691                                switch((int)iomodel->control_type[i]){
    682692                                        case RheologyBbarEnum:
    … …  
    684694                                                        _assert_(iomodel->rheology_B);_assert_(iomodel->cm_min); _assert_(iomodel->cm_max);
    685695                                                        for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->rheology_B[int(iomodel->elements[num_vertices*index+j]-1)];
    686                                                         for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->cm_min[int(iomodel->elements[num_vertices*index+j]-1)*iomodel->num_control_type+i];
    687                                                         for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->cm_max[int(iomodel->elements[num_vertices*index+j]-1)*iomodel->num_control_type+i];
     696                                                        for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->cm_min[int(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
     697                                                        for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->cm_max[int(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
    688698                                                        this->inputs->AddInput(new ControlInput(RheologyBbarEnum,TriaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    689699                                                }
    … …  
    695705
    696706        /*if 3d*/
    697         else if(iomodel->dim==3){
     707        else if(dim==3){
    698708
    699709                /*Intermediaries*/
    … …  
    716726
    717727                /*Control Inputs*/
    718                 if (iomodel->control_analysis && iomodel->control_type){
    719                         for(i=0;i<iomodel->num_control_type;i++){
     728                if (control_analysis && iomodel->control_type){
     729                        for(i=0;i<num_control_type;i++){
    720730                                switch((int)iomodel->control_type[i]){
    721731                                        case RheologyBbarEnum:
    … …  
    723733                                                        _assert_(iomodel->rheology_B);_assert_(iomodel->cm_min); _assert_(iomodel->cm_max);
    724734                                                        for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->rheology_B[int(iomodel->elements[num_vertices*index+j]-1)];
    725                                                         for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->cm_min[int(iomodel->elements[num_vertices*index+j]-1)*iomodel->num_control_type+i];
    726                                                         for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->cm_max[int(iomodel->elements[num_vertices*index+j]-1)*iomodel->num_control_type+i];
     735                                                        for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->cm_min[int(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
     736                                                        for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->cm_max[int(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
    727737                                                        this->inputs->AddInput(new ControlInput(RheologyBEnum,PentaVertexInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    728738                                                }

  • issm/trunk/src/c/objects/Materials/Matpar.cpp

    r9320 r9356  
    2626
    2727        this->mid                       = matpar_mid;
    28         this->rho_ice                   = iomodel->rho_ice;
    29         this->rho_water                 = iomodel->rho_water;
    30         this->heatcapacity              = iomodel->heatcapacity;
    31         this->thermalconductivity       = iomodel->thermalconductivity;
    32         this->latentheat                = iomodel->latentheat;
    33         this->beta                      = iomodel->beta;
    34         this->meltingpoint              = iomodel->meltingpoint;
    35         this->referencetemperature      = iomodel->referencetemperature;
    36         this->mixed_layer_capacity      = iomodel->mixed_layer_capacity;
    37         this->thermal_exchange_velocity = iomodel->thermal_exchange_velocity;
    38         this->g                         = iomodel->g;
    39 
    40         this->kn                        = iomodel->hydro_kn;                     
    41         this->hydro_p                   = iomodel->hydro_p;             
    42         this->hydro_q                   = iomodel->hydro_q;             
    43         this->hydro_CR                  = iomodel->hydro_CR;                     
    44         this->hydro_n                   = iomodel->hydro_n;
     28        iomodel->parameters->FindParam(&this->rho_ice,RhoIceEnum);
     29        iomodel->parameters->FindParam(&this->rho_water,RhoWaterEnum);
     30        iomodel->parameters->FindParam(&this->heatcapacity,HeatcapacityEnum);
     31        iomodel->parameters->FindParam(&this->thermalconductivity,ThermalconductivityEnum);
     32        iomodel->parameters->FindParam(&this->latentheat,LatentheatEnum);
     33        iomodel->parameters->FindParam(&this->beta,BetaEnum);
     34        iomodel->parameters->FindParam(&this->meltingpoint,MeltingpointEnum);
     35        iomodel->parameters->FindParam(&this->referencetemperature,ReferencetemperatureEnum);
     36        iomodel->parameters->FindParam(&this->mixed_layer_capacity,MixedLayerCapacityEnum);
     37        iomodel->parameters->FindParam(&this->thermal_exchange_velocity,ThermalExchangeVelocityEnum);
     38        iomodel->parameters->FindParam(&this->g,GEnum);
     39       
     40        iomodel->parameters->FindParam(&this->kn,HydroKnEnum);
     41        iomodel->parameters->FindParam(&this->hydro_p,HydroPEnum);
     42        iomodel->parameters->FindParam(&this->hydro_q,HydroQEnum);
     43        iomodel->parameters->FindParam(&this->hydro_CR,HydroCREnum);
     44        iomodel->parameters->FindParam(&this->hydro_n,HydroNEnum);
    4545}
    4646/*}}}1*/

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

    r9320 r9356  
    3434        int k,l;
    3535        int gsize;
     36        int dim;
     37
     38        /*Fetch parameters: */
     39        iomodel->parameters->FindParam(&dim,DimEnum);
    3640
    3741        /*id: */
    … …  
    6973        /*Diagnostic Horiz*/
    7074        if (analysis_type==DiagnosticHorizAnalysisEnum){
    71                 if (iomodel->dim==3){
     75                if (dim==3){
    7276                        /*We have a  3d mesh, we may have collapsed elements, hence dead nodes. Freeze them out: */
    7377                        if (!iomodel->nodeonbed) _error_("iomodel->nodeonbed is NULL");
    … …  
    116120                                analysis_type==BalancethicknessAnalysisEnum
    117121                                ){
    118                 if (iomodel->dim==3){
     122                if (dim==3){
    119123                        /*On a 3d mesh, we may have collapsed elements, hence dead nodes. Freeze them out: */
    120124                        _assert_(iomodel->nodeonbed);

  • issm/trunk/src/c/objects/Params/BoolParam.cpp

    r9320 r9356  
    141141#endif
    142142/*}}}*/
     143/*FUNCTION BoolParam::UnitConversion{{{1*/
     144void  BoolParam::UnitConversion(int direction_enum){
     145        /*do nothing, no unit conversion*/
     146}
     147/*}}}*/

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

    r9320 r9356  
    7676                void  SetValue(FILE* fid){_error_("Bool param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7777                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("Bool param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     78                void  UnitConversion(int direction_enum);
    7879               
    7980                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/DoubleMatArrayParam.cpp

    r9320 r9356  
    384384}
    385385/*}}}*/
     386/*FUNCTION DoubleMatArrayParam::UnitConversion{{{1*/
     387void  DoubleMatArrayParam::UnitConversion(int direction_enum){
     388        /*go through all matrices and convert: */
     389        for (int i=0;i<this->M;i++){
     390                double* matrix=this->array[i];
     391                int     m=this->mdim_array[i];
     392                int     n=this->ndim_array[i];
     393                ::UnitConversion(matrix,m*n,direction_enum,this->enum_type);
     394        }
     395
     396}
     397/*}}}*/

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

    r9320 r9356  
    7979                void  SetValue(FILE* fid){_error_("Bool param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    8080                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array);
     81                void  UnitConversion(int direction_enum);
    8182
    8283                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/DoubleMatParam.cpp

    r9320 r9356  
    224224}
    225225/*}}}*/
     226/*FUNCTION DoubleMatParam::UnitConversion{{{1*/
     227void  DoubleMatParam::UnitConversion(int direction_enum){
     228        ::UnitConversion(this->value,this->M*this->N,direction_enum,this->enum_type);
     229}
     230/*}}}*/

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

    r9320 r9356  
    7878                void  SetValue(FILE* fid){_error_("DoubleMat param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7979                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("DoubleMat param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     80                void  UnitConversion(int direction_enum);
    8081
    8182                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/DoubleParam.cpp

    r9320 r9356  
    227227#endif
    228228/*}}}*/
     229/*FUNCTION DoubleParam::UnitConversion{{{1*/
     230void  DoubleParam::UnitConversion(int direction_enum){
     231        ::UnitConversion(&this->value,1,direction_enum,this->enum_type);
     232}
     233/*}}}*/

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

    r9320 r9356  
    7777                void  SetValue(FILE* fid){_error_("Double param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7878                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("Double param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     79                void  UnitConversion(int direction_enum);
    7980
    8081                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/DoubleVecParam.cpp

    r9320 r9356  
    215215}
    216216/*}}}*/
     217/*FUNCTION DoubleVecParam::UnitConversion{{{1*/
     218void  DoubleVecParam::UnitConversion(int direction_enum){
     219        ::UnitConversion(this->values,this->M,direction_enum,this->enum_type);
     220}
     221/*}}}*/

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

    r9320 r9356  
    7777                void  SetValue(FILE* fid){_error_("DoubleVec param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7878                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("DoubleVec param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     79                void  UnitConversion(int direction_enum);
    7980               
    8081                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/FileParam.cpp

    r9320 r9356  
    107107#endif
    108108/*}}}*/
     109/*FUNCTION FileParam::UnitConversion{{{1*/
     110void  FileParam::UnitConversion(int direction_enum){
     111        /*do nothing, no unit conversion*/
     112}
     113/*}}}*/

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

    r9320 r9356  
    7676                void  SetValue(FILE* fid){_error_("File param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7777                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("File param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     78                void  UnitConversion(int direction_enum);
    7879
    7980                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/IntMatParam.cpp

    r9320 r9356  
    212212}
    213213/*}}}*/
     214/*FUNCTION IntMatParam::UnitConversion{{{1*/
     215void  IntMatParam::UnitConversion(int direction_enum){
     216        /*do nothing, no unit conversion*/
     217}
     218/*}}}*/

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

    r9320 r9356  
    7878                void  SetValue(FILE* fid){_error_("IntMat param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7979                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("IntMat param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     80                void  UnitConversion(int direction_enum);
    8081
    8182                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/IntParam.cpp

    r9320 r9356  
    142142#endif
    143143/*}}}*/
     144/*FUNCTION IntParam::UnitConversion{{{1*/
     145void  IntParam::UnitConversion(int direction_enum){
     146        /*do nothing, no unit conversion*/
     147}
     148/*}}}*/

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

    r9320 r9356  
    7777                void  SetValue(FILE* fid){_error_("Int param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7878                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("Int param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     79                void  UnitConversion(int direction_enum);
    7980
    8081                char* GetParameterName(void);

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

    r9320 r9356  
    227227}
    228228/*}}}*/
     229/*FUNCTION IntVecParam::UnitConversion{{{1*/
     230void  IntVecParam::UnitConversion(int direction_enum){
     231        /*do nothing, no unit conversion*/
     232}
     233/*}}}*/

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

    r9320 r9356  
    7878                void  SetValue(FILE* fid){_error_("IntVec param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7979                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("IntVec param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     80                void  UnitConversion(int direction_enum);
    8081               
    8182                char* GetParameterName(void);

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

    r9320 r9356  
    5858                virtual void  SetValue(FILE* fid)=0;
    5959                virtual void  SetValue(double** array, int M, int* mdim_array, int* ndim_array)=0;
     60                virtual void  UnitConversion(int direction_enum)=0;
    6061
    6162                virtual char* GetParameterName(void)=0;

  • issm/trunk/src/c/objects/Params/PetscMatParam.cpp

    r9320 r9356  
    234234}
    235235/*}}}*/
     236/*FUNCTION PetscMatParam::UnitConversion{{{1*/
     237void  PetscMatParam::UnitConversion(int direction_enum){
     238        /*do nothing, no unit conversion*/
     239}
     240/*}}}*/

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

    r9320 r9356  
    7777                void  SetValue(FILE* fid){_error_("PetscMat param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7878                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("PetscMat param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     79                void  UnitConversion(int direction_enum);
    7980
    8081                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/PetscVecParam.cpp

    r9320 r9356  
    229229}
    230230/*}}}*/
     231/*FUNCTION PetscVecParam::UnitConversion{{{1*/
     232void  PetscVecParam::UnitConversion(int direction_enum){
     233        /*do nothing, no unit conversion*/
     234}
     235/*}}}*/

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

    r9320 r9356  
    7777                void  SetValue(FILE* fid){_error_("PetscVec of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7878                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("PetscVec param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     79                void  UnitConversion(int direction_enum);
    7980
    8081                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/StringArrayParam.cpp

    r9320 r9356  
    285285}
    286286/*}}}*/
     287/*FUNCTION StringArrayParam::UnitConversion{{{1*/
     288void  StringArrayParam::UnitConversion(int direction_enum){
     289        /*do nothing, no unit conversion*/
     290}
     291/*}}}*/

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

    r9320 r9356  
    7979                void  SetValue(FILE* fid){_error_("StringArray param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    8080                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("StringArray param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     81                void  UnitConversion(int direction_enum);
    8182
    8283                char* GetParameterName(void);

  • issm/trunk/src/c/objects/Params/StringParam.cpp

    r9336 r9356  
    188188}
    189189/*}}}*/
     190/*FUNCTION StringParam::UnitConversion{{{1*/
     191void  StringParam::UnitConversion(int direction_enum){
     192        /*do nothing, no unit conversion*/
     193}
     194/*}}}*/

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

    r9320 r9356  
    7777                void  SetValue(FILE* fid){_error_("String param of enum %i (%s) cannot hold a FILE",enum_type,EnumToStringx(enum_type));}
    7878                void  SetValue(double** array, int M, int* mdim_array, int* ndim_array){_error_("String param of enum %i (%s) cannot hold an array of matrices",enum_type,EnumToStringx(enum_type));}
     79                void  UnitConversion(int direction_enum);
    7980
    8081                char* GetParameterName(void);

  • issm/trunk/src/c/shared/Numerics/UnitConversion.cpp

    r9320 r9356  
    4747        double scale;
    4848        switch(type_enum){
     49                case DtEnum:        scale=1.0/yts;break; //yr
     50                case NdtEnum:        scale=1.0/yts;break; //yr
    4951                case TimeEnum:        scale=1.0/yts;break; //yr
    5052                case VxEnum:          scale=yts;break; //m/yr

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

    r8926 r9356  
    4646
    4747        /*Recover parameters used throughout the solution:{{{1*/
    48         femmodel->parameters->FindParam(&num_controls,NumControlsEnum);
    49         femmodel->parameters->FindParam(&num_responses,NumResponsesEnum);
     48        femmodel->parameters->FindParam(&num_controls,NumControlTypeEnum);
     49        femmodel->parameters->FindParam(&num_responses,NumCmResponsesEnum);
    5050        femmodel->parameters->FindParam(&control_type,NULL,ControlTypeEnum);
    5151        femmodel->parameters->FindParam(&responses,NULL,NULL,CmResponsesEnum);

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

    r8926 r9356  
    1515
    1616        /*retrieve output file name: */
    17         femmodel->parameters->FindParam(&num_controls,NumControlsEnum);
     17        femmodel->parameters->FindParam(&num_controls,NumControlTypeEnum);
    1818        femmodel->parameters->FindParam(&control_type,NULL,ControlTypeEnum);
    1919        femmodel->parameters->FindParam(&nsteps,NstepsEnum);

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

    r8926 r9356  
    2929        /*retrieve parameters:*/
    3030        femmodel->parameters->FindParam(&control_steady,ControlSteadyEnum);
    31         femmodel->parameters->FindParam(&num_controls,NumControlsEnum);
     31        femmodel->parameters->FindParam(&num_controls,NumControlTypeEnum);
    3232        femmodel->parameters->FindParam(&control_type,NULL,ControlTypeEnum);
    3333        femmodel->parameters->FindParam(&optscal_list,NULL,NULL,OptscalEnum);

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

    r9219 r9356  
    2121        int    solution_type;
    2222        int    output_frequency;
    23         int    dim,gl_migration;
     23        int    dim,groundingline_migration;
    2424       
    2525        /*intermediary: */
    … …  
    3636        femmodel->parameters->FindParam(&output_frequency,OutputFrequencyEnum);
    3737        femmodel->parameters->FindParam(&time_adapt,TimeAdaptEnum);
    38         femmodel->parameters->FindParam(&gl_migration,GroundingLineMigrationEnum);
     38        femmodel->parameters->FindParam(&groundingline_migration,GroundinglineMigrationEnum);
    3939        femmodel->parameters->FindParam(&isdiagnostic,IsdiagnosticEnum);
    4040        femmodel->parameters->FindParam(&isprognostic,IsprognosticEnum);
    … …  
    7676                }
    7777
    78                 if (gl_migration!=NoneEnum){
     78                if (groundingline_migration!=NoneEnum){
    7979                        if(dim==3) _error_("Grounding line migration not implemented in 3d");
    8080                        _printf_(VerboseSolution(),"   computing new grounding line position\n");
Note: See TracChangeset for help on using the changeset viewer.