Changeset 13129


Ignore:
Timestamp:
08/22/12 15:38:04 (13 years ago)
Author:
Mathieu Morlighem
Message:

NEW: tria->material is now an object-oriented that can either be Matice or Matdamageice

Location:
issm/trunk-jpl/src/c
Files:
22 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/EnumDefinitions/EnumDefinitions.h

    r13119 r13129  
    282282        MacAyeal3dIceFrontEnum,
    283283        MaticeEnum,
     284        MatdamageiceEnum,
    284285        MatparEnum,
    285286        NodeEnum,
  • issm/trunk-jpl/src/c/Makefile.am

    r13086 r13129  
    109109                                        ./classes/objects/Materials/Matice.h\
    110110                                        ./classes/objects/Materials/Matice.cpp\
     111                                        ./classes/objects/Materials/Matdamageice.h\
     112                                        ./classes/objects/Materials/Matdamageice.cpp\
    111113                                        ./classes/objects/Materials/Matpar.h\
    112114                                        ./classes/objects/Materials/Matpar.cpp\
  • issm/trunk-jpl/src/c/classes/objects/Elements/Penta.cpp

    r13119 r13129  
    3030
    3131        this->nodes=NULL;
    32         this->matice=NULL;
     32        this->material=NULL;
    3333        this->matpar=NULL;
    3434        this->verticalneighbors=NULL;
     
    4949Penta::Penta(int penta_id, int penta_sid, int index, IoModel* iomodel,int nummodels)
    5050        :PentaRef(nummodels)
    51         ,PentaHook(nummodels,index+1,iomodel) //index+1: matice id, iomodel->numberofelements+1: matpar id
     51        ,PentaHook(nummodels,index+1,iomodel) //index+1: material id, iomodel->numberofelements+1: matpar id
    5252                                                                      { //i is the element index
    5353
     
    8585        /*initialize pointers:*/
    8686        this->nodes=NULL;
    87         this->matice=NULL;
     87        this->material=NULL;
    8888        this->matpar=NULL;
    8989        this->verticalneighbors=NULL;
     
    107107        penta->hnodes=new Hook*[penta->numanalyses];
    108108        for(i=0;i<penta->numanalyses;i++)penta->hnodes[i]=(Hook*)this->hnodes[i]->copy();
    109         penta->hmatice=(Hook*)this->hmatice->copy();
     109        penta->hmaterial=(Hook*)this->hmaterial->copy();
    110110        penta->hmatpar=(Hook*)this->hmatpar->copy();
    111111        penta->hneighbors=(Hook*)this->hneighbors->copy();
     
    132132        penta->nodes=xNew<Node*>(6); //we cannot rely on an analysis_counter to tell us which analysis_type we are running, so we just copy the nodes.
    133133        for(i=0;i<6;i++)penta->nodes[i]=this->nodes[i];
    134         penta->matice=(Matice*)penta->hmatice->delivers();
     134        penta->material=(Material*)penta->hmaterial->delivers();
    135135        penta->matpar=(Matpar*)penta->hmatpar->delivers();
    136136        penta->verticalneighbors=(Penta**)penta->hneighbors->deliverp();
     
    286286                /*Compute strain rate viscosity and pressure: */
    287287                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    288                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     288                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    289289                pressure_input->GetInputValue(&pressure,gauss);
    290290
     
    353353                /*Compute strain rate viscosity and pressure: */
    354354                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    355                 matice->GetViscosity3d(&viscosity,&epsilon[0]);
     355                material->GetViscosity3d(&viscosity,&epsilon[0]);
    356356                pressure_input->GetInputValue(&pressure,gauss);
    357357
     
    391391         * datasets, using internal ids and offsets hidden in hooks: */
    392392        if (this->hnodes[analysis_counter]) this->hnodes[analysis_counter]->configure(nodesin);
    393         this->hmatice->configure(materialsin);
     393        this->hmaterial->configure(materialsin);
    394394        this->hmatpar->configure(materialsin);
    395395        this->hneighbors->configure(elementsin);
     
    398398        if (this->hnodes[analysis_counter]) this->nodes=(Node**)this->hnodes[analysis_counter]->deliverp();
    399399        else this->nodes=NULL;
    400         this->matice=(Matice*)this->hmatice->delivers();
     400        this->material=(Material*)this->hmaterial->delivers();
    401401        this->matpar=(Matpar*)this->hmatpar->delivers();
    402402        this->verticalneighbors=(Penta**)this->hneighbors->deliverp();
     
    419419
    420420        /*Checks in debugging {{{*/
    421         _assert_(this->nodes && this->matice && this->matpar && this->verticalneighbors && this->parameters && this->inputs);
     421        _assert_(this->nodes && this->material && this->matpar && this->verticalneighbors && this->parameters && this->inputs);
    422422        /*}}}*/
    423423       
     
    490490        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    491491        ElementMatrix* Ke=tria->CreateKMatrixPrognostic();
    492         delete tria->matice; delete tria;
     492        delete tria->material; delete tria;
    493493
    494494        /*Delete Vx and Vy averaged*/
     
    507507        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    508508        ElementMatrix* Ke=tria->CreateMassMatrix();
    509         delete tria->matice; delete tria;
     509        delete tria->material; delete tria;
    510510
    511511        /*clean up and return*/
     
    522522
    523523        /*if debugging mode, check that all pointers exist {{{*/
    524         _assert_(this->nodes && this->matice && this->matpar && this->verticalneighbors && this->parameters && this->inputs);
     524        _assert_(this->nodes && this->material && this->matpar && this->verticalneighbors && this->parameters && this->inputs);
    525525        /*}}}*/
    526526
     
    591591        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    592592        ElementVector* pe=tria->CreatePVectorPrognostic();
    593         delete tria->matice; delete tria;
     593        delete tria->material; delete tria;
    594594
    595595        /*Delete Vx and Vy averaged*/
     
    609609        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    610610        ElementVector* pe=tria->CreatePVectorSlope();
    611         delete tria->matice; delete tria;
     611        delete tria->material; delete tria;
    612612
    613613        /*clean up and return*/
     
    624624
    625625        /*Checks in debugging {{{*/
    626         _assert_(this->nodes && this->matice && this->matpar && this->verticalneighbors && this->parameters && this->inputs);
     626        _assert_(this->nodes && this->material && this->matpar && this->verticalneighbors && this->parameters && this->inputs);
    627627        /*}}}*/
    628628
     
    661661        nodes[4]->DeepEcho();
    662662        nodes[5]->DeepEcho();
    663         matice->DeepEcho();
     663        material->DeepEcho();
    664664        matpar->DeepEcho();
    665665        _printLine_("   neighbor ids: " << verticalneighbors[0]->Id() << "-" << verticalneighbors[1]->Id());
     
    13881388                 original_input=(Input*)penta->inputs->GetInput(enum_type);
    13891389                else if (object_enum==MaterialsEnum)
    1390                  original_input=(Input*)penta->matice->inputs->GetInput(enum_type);
     1390                 original_input=(Input*)penta->material->inputs->GetInput(enum_type);
    13911391                else
    13921392                 _error_("object " << EnumToStringx(object_enum) << " not supported yet");
     
    14481448         this->inputs->AddInput((Input*)depth_averaged_input);
    14491449        else if (object_enum==MaterialsEnum)
    1450          this->matice->inputs->AddInput((Input*)depth_averaged_input);
     1450         this->material->inputs->AddInput((Input*)depth_averaged_input);
    14511451        else
    14521452         _error_("object " << EnumToStringx(object_enum) << " not supported yet");
     
    14811481                num_inputs=1;
    14821482                base_inputs=xNew<Input*>(num_inputs);
    1483                 base_inputs[0]=(Input*)matice->inputs->GetInput(enum_type);
     1483                base_inputs[0]=(Input*)material->inputs->GetInput(enum_type);
    14841484        }
    14851485        else if (object_type==NodeEnum){
     
    15151515                        }
    15161516                        else if(object_type==MaterialsEnum){
    1517                                 penta->matice->inputs->AddInput((Input*)copy);
     1517                                penta->material->inputs->AddInput((Input*)copy);
    15181518                        }
    15191519                        else if(object_type==NodeEnum){
     
    15541554
    15551555        /*Go through all the input objects, and find the one corresponding to enum_type, if it exists: */
    1556         if (enum_type==MaterialsRheologyBbarEnum) input=this->matice->inputs->GetInput(MaterialsRheologyBEnum);
    1557         else if (enum_type==MaterialsRheologyZbarEnum) input=this->matice->inputs->GetInput(MaterialsRheologyZEnum);
     1556        if (enum_type==MaterialsRheologyBbarEnum) input=this->material->inputs->GetInput(MaterialsRheologyBEnum);
     1557        else if (enum_type==MaterialsRheologyZbarEnum) input=this->material->inputs->GetInput(MaterialsRheologyZEnum);
    15581558        else input=this->inputs->GetInput(enum_type);
    15591559        //if (!input) _error_("Input " << EnumToStringx(enum_type) << " not found in penta->inputs"); why error out? if the requested input does not exist, we should still
     
    16691669                                case MaterialsRheologyBbarEnum:
    16701670                                case MaterialsRheologyZbarEnum:
    1671                                         /*Matice will take care of it*/ break;
     1671                                        /*Material will take care of it*/ break;
    16721672                                default:
    16731673                                        _error_("Control " << EnumToStringx(reCast<int,IssmDouble>(iomodel->Data(InversionControlParametersEnum)[i])) << " not implemented yet");
     
    19651965                                /*update input*/
    19661966                                if (name==MaterialsRheologyZEnum || name==MaterialsRheologyZbarEnum){
    1967                                         matice->inputs->AddInput(new PentaP1Input(name,values));
     1967                                        material->inputs->AddInput(new PentaP1Input(name,values));
    19681968                                }
    19691969                                else{
     
    26372637        this->SpawnTriaHook(dynamic_cast<TriaHook*>(tria),&indices[0]);
    26382638
    2639         /*Spawn matice*/
    2640         tria->matice=NULL;
    2641         tria->matice=(Matice*)this->matice->copy();
    2642         delete tria->matice->inputs;
    2643         tria->matice->inputs=(Inputs*)this->matice->inputs->SpawnTriaInputs(indices);
    2644 
    2645         /*recover nodes, matice and matpar: */
     2639        /*Spawn material*/
     2640        tria->material=NULL;
     2641        tria->material=(Material*)this->material->copy();
     2642        delete tria->material->inputs;
     2643        tria->material->inputs=(Inputs*)this->material->inputs->SpawnTriaInputs(indices);
     2644
     2645        /*recover nodes, material and matpar: */
    26462646        tria->nodes=(Node**)tria->hnodes[analysis_counter]->deliverp();
    26472647        tria->matpar=(Matpar*)tria->hmatpar->delivers();
     
    27302730                tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
    27312731                S=tria->SurfaceArea();
    2732                 delete tria->matice; delete tria;
     2732                delete tria->material; delete tria;
    27332733                return S;
    27342734        }
     
    27372737                tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
    27382738                S=tria->SurfaceArea();
    2739                 delete tria->matice; delete tria;
     2739                delete tria->material; delete tria;
    27402740                return S;
    27412741        }
     
    30183018
    30193019                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    3020                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     3020                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    30213021                GetPhi(&phi, &epsilon[0], viscosity);
    30223022               
     
    31283128        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    31293129        mass_flux=tria->MassFlux(segment,process_units);
    3130         delete tria->matice; delete tria;
     3130        delete tria->material; delete tria;
    31313131
    31323132        /*Delete Vx and Vy averaged*/
     
    32353235        switch(response_enum){
    32363236                case MaterialsRheologyBbarEnum:
    3237                         *presponse=this->matice->GetBbar();
     3237                        *presponse=this->material->GetBbar();
    32383238                        break;
    32393239                case MaterialsRheologyZbarEnum:
    3240                         *presponse=this->matice->GetZbar();
     3240                        *presponse=this->material->GetZbar();
    32413241                        break;
    32423242                case VelEnum:
     
    35313531        ElementMatrix* Ke=tria->CreateKMatrixMelting();
    35323532
    3533         delete tria->matice; delete tria;
     3533        delete tria->material; delete tria;
    35343534        return Ke;
    35353535}
     
    38303830
    38313831                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    3832                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     3832                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    38333833                GetPhi(&phi, &epsilon[0], viscosity);
    38343834
     
    40864086
    40874087                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    4088                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     4088                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    40894089                GetPhi(&phi, &epsilon[0], viscosity);
    40904090
     
    43564356                                B_average=Paterson((values[0]+values[1]+values[2]+values[3]+values[4]+values[5])/6.0);
    43574357                                for(i=0;i<numdof;i++) B[i]=B_average;
    4358                                 this->matice->inputs->AddInput(new PentaP1Input(MaterialsRheologyBEnum,B));
     4358                                this->material->inputs->AddInput(new PentaP1Input(MaterialsRheologyBEnum,B));
    43594359                                break;
    43604360                        case ArrheniusEnum:
     
    43624362                                B_average=Arrhenius((values[0]+values[1]+values[2]+values[3]+values[4]+values[5])/6.0,
    43634363                                                        s_average-((xyz_list[0][2]+xyz_list[1][2]+xyz_list[2][2]+xyz_list[3][2]+xyz_list[4][2]+xyz_list[5][2])/6.0),
    4364                                                         matice->GetN());
     4364                                                        material->GetN());
    43654365                                for(i=0;i<numdof;i++) B[i]=B_average;
    4366                                 this->matice->inputs->AddInput(new PentaP1Input(MaterialsRheologyBEnum,B));
     4366                                this->material->inputs->AddInput(new PentaP1Input(MaterialsRheologyBEnum,B));
    43674367                                break;
    43684368                        default:
     
    44344434                                B_average=Paterson((temperatures[0]+temperatures[1]+temperatures[2]+temperatures[3]+temperatures[4]+temperatures[5])/6.0);
    44354435                                for(i=0;i<numdof;i++) B[i]=B_average;
    4436                                 this->matice->inputs->AddInput(new PentaP1Input(MaterialsRheologyBEnum,B));
     4436                                this->material->inputs->AddInput(new PentaP1Input(MaterialsRheologyBEnum,B));
    44374437                                break;
    44384438                        case ArrheniusEnum:
     
    44404440                                B_average=Arrhenius((temperatures[0]+temperatures[1]+temperatures[2]+temperatures[3]+temperatures[4]+temperatures[5])/6.0,
    44414441                                                        s_average-((xyz_list[0][2]+xyz_list[1][2]+xyz_list[2][2]+xyz_list[3][2]+xyz_list[4][2]+xyz_list[5][2])/6.0),
    4442                                                         matice->GetN());
     4442                                                        material->GetN());
    44434443                                for(i=0;i<numdof;i++) B[i]=B_average;
    4444                                 this->matice->inputs->AddInput(new PentaP1Input(MaterialsRheologyBEnum,B));
     4444                                this->material->inputs->AddInput(new PentaP1Input(MaterialsRheologyBEnum,B));
    44454445                                break;
    44464446                        default:
     
    44684468        if(enum_type==MaterialsRheologyBbarEnum){
    44694469                if(!IsOnBed()) return;
    4470                 input=(Input*)matice->inputs->GetInput(MaterialsRheologyBEnum);
     4470                input=(Input*)material->inputs->GetInput(MaterialsRheologyBEnum);
    44714471        }
    44724472        else if(enum_type==MaterialsRheologyZbarEnum){
    44734473                if(!IsOnBed()) return;
    4474                 input=(Input*)matice->inputs->GetInput(MaterialsRheologyZEnum);
     4474                input=(Input*)material->inputs->GetInput(MaterialsRheologyZEnum);
    44754475        }
    44764476               
     
    44914491
    44924492        if(enum_type==MaterialsRheologyBbarEnum){
    4493                 input=(Input*)matice->inputs->GetInput(MaterialsRheologyBEnum);
     4493                input=(Input*)material->inputs->GetInput(MaterialsRheologyBEnum);
    44944494        }
    44954495        else if(enum_type==MaterialsRheologyZbarEnum){
    4496                 input=(Input*)matice->inputs->GetInput(MaterialsRheologyZEnum);
     4496                input=(Input*)material->inputs->GetInput(MaterialsRheologyZEnum);
    44974497        }
    44984498        else{
     
    45134513
    45144514        if(enum_type==MaterialsRheologyBbarEnum){
    4515                 input=(Input*)matice->inputs->GetInput(MaterialsRheologyBEnum);
     4515                input=(Input*)material->inputs->GetInput(MaterialsRheologyBEnum);
    45164516        }
    45174517        else if(enum_type==MaterialsRheologyZbarEnum){
    4518                 input=(Input*)matice->inputs->GetInput(MaterialsRheologyZEnum);
     4518                input=(Input*)material->inputs->GetInput(MaterialsRheologyZEnum);
    45194519        }
    45204520        else{
     
    45584558        if (!IsOnBed()) return NULL;
    45594559
    4560         /*Depth Averaging B*/
    4561         this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    4562 
    4563         /*Depth Averaging Z*/
    4564         this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     4560        /*Depth average some fields*/
     4561        switch(this->material->ObjectEnum()){
     4562                case MaticeEnum:
     4563                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     4564                        break;
     4565                case MatdamageiceEnum:
     4566                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     4567                        this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     4568                        break;
     4569                default:
     4570                        _error_("material "<<EnumToStringx(this->material->ObjectEnum())<<" not supported");
     4571        }
    45654572
    45664573        /*Call Tria function*/
    45674574        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    45684575        ElementMatrix* Ke=tria->CreateKMatrixAdjointMacAyeal();
    4569         delete tria->matice; delete tria;
    4570 
    4571         /*Delete B averaged*/
    4572         this->matice->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    4573 
    4574         /*Delete Z averaged*/
    4575         this->matice->inputs->DeleteInput(MaterialsRheologyZbarEnum);
     4576        delete tria->material; delete tria;
     4577
     4578        /*Delete averaged fields*/
     4579        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     4580        this->material->inputs->DeleteInput(MaterialsRheologyZbarEnum);
    45764581
    45774582        /*clean up and return*/
     
    46194624
    46204625                this->GetStrainRate3dPattyn(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    4621                 matice->GetViscosityDerivativeEpsSquare(&mu_prime,&epsilon[0]);
     4626                material->GetViscosityDerivativeEpsSquare(&mu_prime,&epsilon[0]);
    46224627                eps1[0]=2*epsilon[0]+epsilon[1];   eps2[0]=epsilon[2];
    46234628                eps1[1]=epsilon[2];                eps2[1]=epsilon[0]+2*epsilon[1];
     
    46894694
    46904695                this->GetStrainRate3dPattyn(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    4691                 matice->GetViscosityDerivativeEpsSquare(&mu_prime,&epsilon[0]);
     4696                material->GetViscosityDerivativeEpsSquare(&mu_prime,&epsilon[0]);
    46924697                eps1[0]=epsilon[0];   eps2[0]=epsilon[2];   eps3[0]=epsilon[3];
    46934698                eps1[1]=epsilon[2];   eps2[1]=epsilon[1];   eps3[1]=epsilon[4];
     
    47544759        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    47554760        ElementVector* pe=tria->CreatePVectorAdjointHoriz();
    4756         delete tria->matice; delete tria;
     4761        delete tria->material; delete tria;
    47574762
    47584763        /*clean up and return*/
     
    47684773        Tria* tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
    47694774        ElementVector* pe=tria->CreatePVectorAdjointHoriz();
    4770         delete tria->matice; delete tria;
     4775        delete tria->material; delete tria;
    47714776
    47724777        /*clean up and return*/
     
    47824787        Tria* tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
    47834788        ElementVector* pe=tria->CreatePVectorAdjointStokes();
    4784         delete tria->matice; delete tria;
     4789        delete tria->material; delete tria;
    47854790
    47864791        /*clean up and return*/
     
    48794884                        tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    48804885                        tria->GradjDragGradient(gradient,resp,control_index);
    4881                         delete tria->matice; delete tria;
     4886                        delete tria->material; delete tria;
    48824887                        break;
    48834888                case RheologyBbarAbsGradientEnum:
     
    48854890                        tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    48864891                        tria->GradjBGradient(gradient,resp,control_index);
    4887                         delete tria->matice; delete tria;
     4892                        delete tria->material; delete tria;
    48884893                        break;
    48894894                default:
     
    49024907        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    49034908        tria->GradjDragMacAyeal(gradient,control_index);
    4904         delete tria->matice; delete tria;
     4909        delete tria->material; delete tria;
    49054910
    49064911} /*}}}*/
     
    50805085        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
    50815086        tria->GradjBMacAyeal(gradient,control_index);
    5082         delete tria->matice; delete tria;
     5087        delete tria->material; delete tria;
    50835088
    50845089        /*delete Average B*/
    5085         this->matice->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     5090        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    50865091
    50875092} /*}}}*/
     
    50985103        Tria* tria=(Tria*)SpawnTria(0,1,2);
    50995104        tria->GradjBMacAyeal(gradient,control_index);    //We use MacAyeal as an estimate for now
    5100         delete tria->matice; delete tria;
     5105        delete tria->material; delete tria;
    51015106
    51025107        /*delete Average B*/
    5103         this->matice->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     5108        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    51045109} /*}}}*/
    51055110/*FUNCTION Penta::GradjBbarStokes {{{*/
     
    51155120        Tria* tria=(Tria*)SpawnTria(0,1,2);
    51165121        tria->GradjBMacAyeal(gradient,control_index);    //We use MacAyeal as an estimate for now
    5117         delete tria->matice; delete tria;
     5122        delete tria->material; delete tria;
    51185123
    51195124        /*delete Average B*/
    5120         this->matice->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     5125        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    51215126} /*}}}*/
    51225127/*FUNCTION Penta::InputControlUpdate{{{*/
     
    51365141                if(control_type[i]==MaterialsRheologyBbarEnum){
    51375142                        if (!IsOnBed()) goto cleanup_and_return;
    5138                         input=(Input*)matice->inputs->GetInput(MaterialsRheologyBEnum); _assert_(input);
     5143                        input=(Input*)material->inputs->GetInput(MaterialsRheologyBEnum); _assert_(input);
    51395144                }
    51405145                else if(control_type[i]==MaterialsRheologyZbarEnum){
    51415146                        if (!IsOnBed()) goto cleanup_and_return;
    5142                         input=(Input*)matice->inputs->GetInput(MaterialsRheologyZEnum); _assert_(input);
     5147                        input=(Input*)material->inputs->GetInput(MaterialsRheologyZEnum); _assert_(input);
    51435148                }
    51445149                else{
     
    52685273                tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
    52695274                J=tria->SurfaceAverageVelMisfit(process_units,weight_index);
    5270                 delete tria->matice; delete tria;
     5275                delete tria->material; delete tria;
    52715276                return J;
    52725277        }
     
    52755280                tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
    52765281                J=tria->SurfaceAverageVelMisfit(process_units,weight_index);
    5277                 delete tria->matice; delete tria;
     5282                delete tria->material; delete tria;
    52785283                return J;
    52795284        }
     
    53055310                tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
    53065311                J=tria->SurfaceAbsVelMisfit(process_units,weight_index);
    5307                 delete tria->matice; delete tria;
     5312                delete tria->material; delete tria;
    53085313                return J;
    53095314        }
     
    53125317                tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
    53135318                J=tria->SurfaceAbsVelMisfit(process_units,weight_index);
    5314                 delete tria->matice; delete tria;
     5319                delete tria->material; delete tria;
    53155320                return J;
    53165321        }
     
    53425347                tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
    53435348                J=tria->SurfaceLogVelMisfit(process_units,weight_index);
    5344                 delete tria->matice; delete tria;
     5349                delete tria->material; delete tria;
    53455350                return J;
    53465351        }
     
    53495354                tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
    53505355                J=tria->SurfaceLogVelMisfit(process_units,weight_index);
    5351                 delete tria->matice; delete tria;
     5356                delete tria->material; delete tria;
    53525357                return J;
    53535358        }
     
    53815386                tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
    53825387                J=tria->SurfaceLogVxVyMisfit(process_units,weight_index);
    5383                 delete tria->matice; delete tria;
     5388                delete tria->material; delete tria;
    53845389                return J;
    53855390        }
     
    53885393                tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
    53895394                J=tria->SurfaceLogVxVyMisfit(process_units,weight_index);
    5390                 delete tria->matice; delete tria;
     5395                delete tria->material; delete tria;
    53915396                return J;
    53925397        }
     
    54185423                tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria (lower face).
    54195424                J=tria->SurfaceRelVelMisfit(process_units,weight_index);
    5420                 delete tria->matice; delete tria;
     5425                delete tria->material; delete tria;
    54215426                return J;
    54225427        }
     
    54255430                tria=(Tria*)SpawnTria(3,4,5); //nodes 3, 4 and 5 make the new tria (upper face).
    54265431                J=tria->SurfaceRelVelMisfit(process_units,weight_index);
    5427                 delete tria->matice; delete tria;
     5432                delete tria->material; delete tria;
    54285433                return J;
    54295434        }
     
    54525457        tria=(Tria*)SpawnTria(0,1,2);
    54535458        J=tria->ThicknessAbsMisfit(process_units,weight_index);
    5454         delete tria->matice; delete tria;
     5459        delete tria->material; delete tria;
    54555460        return J;
    54565461}
     
    54675472        tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria
    54685473        J=tria->DragCoefficientAbsGradient(process_units,weight_index);
    5469         delete tria->matice; delete tria;
     5474        delete tria->material; delete tria;
    54705475        return J;
    54715476}
     
    54825487        tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria
    54835488        J=tria->RheologyBbarAbsGradient(process_units,weight_index);
    5484         delete tria->matice; delete tria;
     5489        delete tria->material; delete tria;
    54855490        return J;
    54865491}
     
    55315536
    55325537        if(control_enum==MaterialsRheologyBbarEnum){
    5533                 input=(Input*)matice->inputs->GetInput(control_enum); _assert_(input);
     5538                input=(Input*)material->inputs->GetInput(control_enum); _assert_(input);
    55345539        }
    55355540        else{
     
    58265831                this->GetStrainRate3dPattyn(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    58275832                this->GetStrainRate3dPattyn(&oldepsilon[0],&xyz_list[0][0],gauss,vxold_input,vyold_input);
    5828                 matice->GetViscosity3d(&viscosity, &epsilon[0]);
    5829                 matice->GetViscosity3d(&oldviscosity, &oldepsilon[0]);
     5833                material->GetViscosity3d(&viscosity, &epsilon[0]);
     5834                material->GetViscosity3d(&oldviscosity, &oldepsilon[0]);
    58305835
    58315836                newviscosity=viscosity+viscosity_overshoot*(viscosity-oldviscosity);
     
    58475852
    58485853        /*Clean-up and return*/
    5849         delete tria->matice; delete tria;
     5854        delete tria->material; delete tria;
    58505855        delete gauss;
    58515856        delete gauss_tria;
     
    60366041
    60376042                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    6038                 matice->GetViscosity3dStokes(&viscosity, &epsilon[0]);
     6043                material->GetViscosity3dStokes(&viscosity, &epsilon[0]);
    60396044
    60406045                D_scalar=2*viscosity*gauss->weight*Jdet;
     
    60636068
    60646069        /*Clean-up and return*/
    6065         delete tria->matice; delete tria;
     6070        delete tria->material; delete tria;
    60666071        delete gauss;
    60676072        delete gauss_tria;
     
    61436148
    61446149                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    6145                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     6150                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    61466151
    61476152                BedNormal(&bed_normal[0],xyz_list_tria);
     
    63356340        if (!IsOnBed()) return NULL;
    63366341
    6337         /*Depth Averaging B*/
    6338         this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    6339 
    6340         /*Depth Averaging Z*/
    6341         this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     6342        /*Depth average some fields*/
     6343        switch(this->material->ObjectEnum()){
     6344                case MaticeEnum:
     6345                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     6346                        break;
     6347                case MatdamageiceEnum:
     6348                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     6349                        this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     6350                        break;
     6351                default:
     6352                        _error_("material "<<EnumToStringx(this->material->ObjectEnum())<<" not supported");
     6353        }
    63426354
    63436355        /*Call Tria function*/
    63446356        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    63456357        ElementMatrix* Ke=tria->CreateKMatrixDiagnosticMacAyeal();
    6346         delete tria->matice; delete tria;
    6347 
    6348         /*Delete B averaged*/
    6349         this->matice->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    6350 
    6351         /*Delete Z averaged*/
    6352         this->matice->inputs->DeleteInput(MaterialsRheologyZbarEnum);
     6358        delete tria->material; delete tria;
     6359
     6360        /*Delete averaged fields*/
     6361        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     6362        this->material->inputs->DeleteInput(MaterialsRheologyZbarEnum);
    63536363
    63546364        /*clean up and return*/
     
    64256435                        this->GetStrainRate3dPattyn(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    64266436                        this->GetStrainRate3dPattyn(&oldepsilon[0],&xyz_list[0][0],gauss,vxold_input,vyold_input);
    6427                         matice->GetViscosity3d(&viscosity, &epsilon[0]);
    6428                         matice->GetViscosity3d(&oldviscosity, &oldepsilon[0]);
     6437                        material->GetViscosity3d(&viscosity, &epsilon[0]);
     6438                        material->GetViscosity3d(&oldviscosity, &oldepsilon[0]);
    64296439
    64306440                        newviscosity=viscosity+viscosity_overshoot*(viscosity-oldviscosity);
     
    64326442                else if (approximation==MacAyealStokesApproximationEnum){
    64336443                        this->GetStrainRate3d(&epsilons[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    6434                         matice->GetViscosity3dStokes(&newviscosity,&epsilons[0]);
     6444                        material->GetViscosity3dStokes(&newviscosity,&epsilons[0]);
    64356445                }
    64366446                else _error_("approximation " << approximation << " (" << EnumToStringx(approximation) << ") not supported yet");
     
    64516461
    64526462        /*Clean up and return*/
    6453         delete tria->matice;
     6463        delete tria->material;
    64546464        delete tria;
    64556465        delete gauss_tria;
     
    64696479        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    64706480        ElementMatrix* Ke=tria->CreateKMatrixDiagnosticMacAyealFriction();
    6471         delete tria->matice; delete tria;
     6481        delete tria->material; delete tria;
    64726482
    64736483        /*clean-up and return*/
     
    65826592
    65836593        /*Clean up and return*/
    6584         delete tria->matice;
     6594        delete tria->material;
    65856595        delete tria;
    65866596        delete gauss_tria;
     
    66006610        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    66016611        ElementMatrix* Ke=tria->CreateKMatrixDiagnosticMacAyealFriction();
    6602         delete tria->matice; delete tria;
     6612        delete tria->material; delete tria;
    66036613
    66046614        /*clean-up and return*/
     
    66656675                this->GetStrainRate3dPattyn(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    66666676                this->GetStrainRate3dPattyn(&oldepsilon[0],&xyz_list[0][0],gauss,vxold_input,vyold_input);
    6667                 matice->GetViscosity3d(&viscosity, &epsilon[0]);
    6668                 matice->GetViscosity3d(&oldviscosity, &oldepsilon[0]);
     6677                material->GetViscosity3d(&viscosity, &epsilon[0]);
     6678                material->GetViscosity3d(&oldviscosity, &oldepsilon[0]);
    66696679                newviscosity=viscosity+viscosity_overshoot*(viscosity-oldviscosity);
    66706680
     
    68296839
    68306840                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    6831                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     6841                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    68326842
    68336843                D_scalar=gauss->weight*Jdet;
     
    69006910
    69016911                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    6902                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     6912                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    69036913
    69046914                friction->GetAlpha2(&alpha2, gauss,VxEnum,VyEnum,VzEnum);
     
    70887098
    70897099                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    7090                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     7100                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    70917101
    70927102                for(i=0;i<NUMVERTICES;i++){
     
    71627172                BedNormal(&bed_normal[0],xyz_list_tria);
    71637173                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    7164                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     7174                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    71657175                friction->GetAlpha2(&alpha2_gauss, gauss,VxEnum,VyEnum,VzEnum);
    71667176
     
    72397249
    72407250                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    7241                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     7251                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    72427252
    72437253                for(i=0;i<NUMVERTICES;i++){
     
    73137323                BedNormal(&bed_normal[0],xyz_list_tria);
    73147324                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    7315                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     7325                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    73167326                friction->GetAlpha2(&alpha2_gauss, gauss,VxEnum,VyEnum,VzEnum);
    73177327
     
    74357445        rho_ice=matpar->GetRhoIce();
    74367446        gravity=matpar->GetG();
    7437         n=matice->GetN();
    7438         B=matice->GetB();
     7447        n=material->GetN();
     7448        B=material->GetB();
    74397449        Input* thickness_input=inputs->GetInput(ThicknessEnum);  _assert_(thickness_input);
    74407450        Input* surface_input=inputs->GetInput(SurfaceEnum);      _assert_(surface_input);
     
    75047514        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    75057515        ElementVector* pe=tria->CreatePVectorDiagnosticMacAyeal();
    7506         delete tria->matice; delete tria;
     7516        delete tria->material; delete tria;
    75077517
    75087518        /*Clean up and return*/
     
    75187528        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    75197529        ElementVector* pe=tria->CreatePVectorDiagnosticMacAyeal();
    7520         delete tria->matice; delete tria;
     7530        delete tria->material; delete tria;
    75217531
    75227532        /*Clean up and return*/
     
    76367646
    76377647                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
    7638                 matice->GetViscosity3dStokes(&viscosity,&epsilon[0]);
     7648                material->GetViscosity3dStokes(&viscosity,&epsilon[0]);
    76397649
    76407650                for(i=0;i<NUMVERTICES+1;i++){
     
    78977907        if (!IsOnBed()) return NULL;
    78987908
    7899         /*Depth Averaging B*/
    7900         this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    7901 
    7902         /*Depth Averaging Z*/
    7903         this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     7909        /*Depth average some fields*/
     7910        switch(this->material->ObjectEnum()){
     7911                case MaticeEnum:
     7912                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     7913                        break;
     7914                case MatdamageiceEnum:
     7915                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     7916                        this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     7917                        break;
     7918                default:
     7919                        _error_("material "<<EnumToStringx(this->material->ObjectEnum())<<" not supported");
     7920        }
    79047921
    79057922        /*Call Tria function*/
    79067923        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    79077924        ElementMatrix* Ke=tria->CreateJacobianDiagnosticMacayeal();
    7908         delete tria->matice; delete tria;
    7909 
    7910         /*Delete B averaged*/
    7911         this->matice->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    7912 
    7913         /*Delete Z averaged*/
    7914         this->matice->inputs->DeleteInput(MaterialsRheologyZbarEnum);
     7925        delete tria->material; delete tria;
     7926
     7927        /*Delete averaged inputs*/
     7928        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     7929        this->material->inputs->DeleteInput(MaterialsRheologyZbarEnum);
    79157930
    79167931        /*clean up and return*/
     
    79557970
    79567971                this->GetStrainRate3dPattyn(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    7957                 matice->GetViscosityDerivativeEpsSquare(&mu_prime,&epsilon[0]);
     7972                material->GetViscosityDerivativeEpsSquare(&mu_prime,&epsilon[0]);
    79587973                eps1[0]=2*epsilon[0]+epsilon[1];   eps2[0]=epsilon[2];
    79597974                eps1[1]=epsilon[2];                eps2[1]=epsilon[0]+2*epsilon[1];
     
    80228037
    80238038                this->GetStrainRate3dPattyn(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    8024                 matice->GetViscosityDerivativeEpsSquare(&mu_prime,&epsilon[0]);
     8039                material->GetViscosityDerivativeEpsSquare(&mu_prime,&epsilon[0]);
    80258040                eps1[0]=epsilon[0];   eps2[0]=epsilon[2];   eps3[0]=epsilon[3];
    80268041                eps1[1]=epsilon[2];   eps2[1]=epsilon[1];   eps3[1]=epsilon[4];
     
    82648279
    82658280        /*Get A*/
    8266         _assert_(matice->GetN()==3.0);
    8267         A=matice->GetA();
     8281        _assert_(material->GetN()==3.0);
     8282        A=material->GetA();
    82688283
    82698284        /*Solve for tau_perp (http://fr.wikipedia.org/wiki/Méthode_de_Cardan)*/
     
    90859100        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    90869101        ElementMatrix* Ke=tria->CreateKMatrixBalancethickness();
    9087         delete tria->matice; delete tria;
     9102        delete tria->material; delete tria;
    90889103
    90899104        /*Delete Vx and Vy averaged*/
     
    91079122        Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
    91089123        ElementVector* pe=tria->CreatePVectorBalancethickness();
    9109         delete tria->matice; delete tria;
     9124        delete tria->material; delete tria;
    91109125
    91119126        /*Delete Vx and Vy averaged*/
  • issm/trunk-jpl/src/c/classes/objects/Elements/Penta.h

    r13092 r13129  
    1616class IoModel;
    1717class Node;
    18 class Matice;
     18class Material;
    1919class Matpar;
    2020class Tria;
     
    3535
    3636                Node       **nodes;        // 6 nodes
    37                 Matice      *matice;       // 1 material ice
     37                Material    *material;       // 1 material ice
    3838                Matpar      *matpar;       // 1 material parameter
    3939                Penta      **verticalneighbors;   // 2 neighbors: first one under, second one above
  • issm/trunk-jpl/src/c/classes/objects/Elements/PentaHook.cpp

    r12832 r13129  
    2525        numanalyses=UNDEF;
    2626        this->hnodes=NULL;
    27         this->hmatice=NULL;
     27        this->hmaterial=NULL;
    2828        this->hmatpar=NULL;
    2929        this->hneighbors=NULL;
     
    3939        }
    4040        delete [] this->hnodes;
    41         delete hmatice;
     41        delete hmaterial;
    4242        delete hmatpar;
    4343        delete hneighbors;
    4444}
    4545/*}}}*/
    46 /*FUNCTION PentaHook::PentaHook(int in_numanalyses,int matice_id, int matpar_id){{{*/
    47 PentaHook::PentaHook(int in_numanalyses,int matice_id, IoModel* iomodel){
     46/*FUNCTION PentaHook::PentaHook(int in_numanalyses,int material_id, int matpar_id){{{*/
     47PentaHook::PentaHook(int in_numanalyses,int material_id, IoModel* iomodel){
    4848
    4949        /*intermediary: */
     
    5555        this->numanalyses=in_numanalyses;
    5656        this->hnodes=new Hook*[in_numanalyses];
    57         this->hmatice=new Hook(&matice_id,1);
     57        this->hmaterial=new Hook(&material_id,1);
    5858        this->hmatpar=new Hook(&matpar_id,1);
    5959        this->hneighbors=NULL;
     
    9898                }
    9999        }
    100         // do not spawn hmatice. matice will be taken care of by Penta
    101         triahook->hmatice=NULL;
     100        // do not spawn hmaterial. material will be taken care of by Penta
     101        triahook->hmaterial=NULL;
    102102        triahook->hmatpar=(Hook*)this->hmatpar->copy();
    103103}
  • issm/trunk-jpl/src/c/classes/objects/Elements/PentaHook.h

    r12365 r13129  
    1515                int   numanalyses; //number of analysis types
    1616                Hook** hnodes; // 6 nodes for each analysis type
    17                 Hook*  hmatice; // 1 ice material
     17                Hook*  hmaterial; // 1 ice material
    1818                Hook*  hmatpar; // 1 material parameter
    1919                Hook*  hneighbors; // 2 elements, first down, second up
     
    2121                /*FUNCTION constructors, destructors {{{*/
    2222                PentaHook();
    23                 PentaHook(int in_numanalyses,int matice_id, IoModel* iomodel);
     23                PentaHook(int in_numanalyses,int material_id, IoModel* iomodel);
    2424                ~PentaHook();
    2525                void SetHookNodes(int* node_ids,int analysis_counter);
  • issm/trunk-jpl/src/c/classes/objects/Elements/Tria.cpp

    r13119 r13129  
    2929
    3030        this->nodes=NULL;
    31         this->matice=NULL;
     31        this->material=NULL;
    3232        this->matpar=NULL;
    3333        for(i=0;i<3;i++)this->horizontalneighborsids[i]=UNDEF;
     
    6161                /*initialize pointers:*/
    6262                this->nodes=NULL;
    63                 this->matice=NULL;
     63                this->material=NULL;
    6464                this->matpar=NULL;
    6565
     
    8989        tria->hnodes=new Hook*[tria->numanalyses];
    9090        for(i=0;i<tria->numanalyses;i++)tria->hnodes[i]=(Hook*)this->hnodes[i]->copy();
    91         tria->hmatice=(Hook*)this->hmatice->copy();
     91        tria->hmaterial=(Hook*)this->hmaterial->copy();
    9292        tria->hmatpar=(Hook*)this->hmatpar->copy();
    9393
     
    113113        tria->nodes=xNew<Node*>(3); //we cannot rely on an analysis_counter to tell us which analysis_type we are running, so we just copy the nodes.
    114114        for(i=0;i<3;i++)tria->nodes[i]=this->nodes[i];
    115         tria->matice=(Matice*)tria->hmatice->delivers();
     115        tria->material=(Material*)tria->hmaterial->delivers();
    116116        tria->matpar=(Matpar*)tria->hmatpar->delivers();
    117117
     
    173173
    174174        /*Checks in debugging mode{{{*/
    175         _assert_(this->nodes && this->matice && this->matpar && this->parameters && this->inputs);
     175        _assert_(this->nodes && this->material && this->matpar && this->parameters && this->inputs);
    176176        /*}}}*/
    177177       
     
    526526        /*asserts: {{{*/
    527527        /*if debugging mode, check that all pointers exist*/
    528         _assert_(this->nodes && this->matice && this->matpar && this->parameters && this->inputs);
     528        _assert_(this->nodes && this->material && this->matpar && this->parameters && this->inputs);
    529529        /*}}}*/
    530530       
     
    748748
    749749        /*Checks in debugging {{{*/
    750         _assert_(this->nodes && this->matice && this->matpar && this->parameters && this->inputs);
     750        _assert_(this->nodes && this->material && this->matpar && this->parameters && this->inputs);
    751751        /*}}}*/
    752752
     
    812812                /*Compute strain rate viscosity and pressure: */
    813813                this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    814                 matice->GetViscosity2d(&viscosity,&epsilon[0]);
     814                material->GetViscosity2d(&viscosity,&epsilon[0]);
    815815                pressure_input->GetInputValue(&pressure,gauss);
    816816
     
    846846         * datasets, using internal ids and offsets hidden in hooks: */
    847847        if(this->hnodes[analysis_counter]) this->hnodes[analysis_counter]->configure(nodesin);
    848         this->hmatice->configure(materialsin);
     848        this->hmaterial->configure(materialsin);
    849849        this->hmatpar->configure(materialsin);
    850850
     
    852852        if(this->hnodes[analysis_counter]) this->nodes=(Node**)this->hnodes[analysis_counter]->deliverp();
    853853        else this->nodes=NULL;
    854         this->matice=(Matice*)this->hmatice->delivers();
     854        this->material=(Material*)this->hmaterial->delivers();
    855855        this->matpar=(Matpar*)this->hmatpar->delivers();
    856856
     
    875875        else _printLine_("nodes = NULL");
    876876
    877         if (matice) matice->DeepEcho();
    878         else _printLine_("matice = NULL");
     877        if (material) material->DeepEcho();
     878        else _printLine_("material = NULL");
    879879
    880880        if (matpar) matpar->DeepEcho();
     
    987987        else _printLine_("nodes = NULL");
    988988
    989         if (matice) matice->Echo();
    990         else _printLine_("matice = NULL");
     989        if (material) material->Echo();
     990        else _printLine_("material = NULL");
    991991
    992992        if (matpar) matpar->Echo();
     
    13501350         oldinput=(Input*)this->inputs->GetInput(enum_type);
    13511351        else if (object_enum==MaterialsEnum)
    1352          oldinput=(Input*)this->matice->inputs->GetInput(enum_type);
     1352         oldinput=(Input*)this->material->inputs->GetInput(enum_type);
    13531353        else
    13541354         _error_("object " << EnumToStringx(object_enum) << " not supported yet");
     
    13631363         this->inputs->AddInput((Input*)newinput);
    13641364        else if (object_enum==MaterialsEnum)
    1365          this->matice->inputs->AddInput((Input*)newinput);
     1365         this->material->inputs->AddInput((Input*)newinput);
    13661366        else
    13671367         _error_("object " << EnumToStringx(object_enum) << " not supported yet");
     
    13961396
    13971397        /*Go through all the input objects, and find the one corresponding to enum_type, if it exists: */
    1398         if (enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum) input=this->matice->inputs->GetInput(enum_type);
     1398        if (enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum) input=this->material->inputs->GetInput(enum_type);
    13991399        else input=this->inputs->GetInput(enum_type);
    14001400        //if (!input) _error_("Input " << EnumToStringx(enum_type) << " not found in tria->inputs");
     
    15031503                                case MaterialsRheologyBbarEnum:
    15041504                                case MaterialsRheologyZbarEnum:
    1505                                         /*Matice will take care of it*/ break;
     1505                                        /*Material will take care of it*/ break;
    15061506                                default:
    15071507                                        _error_("Control " << EnumToStringx(reCast<int,IssmDouble>(iomodel->Data(InversionControlParametersEnum)[i])) << " not implemented yet");
     
    16921692                        /*update input*/
    16931693                        if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheologyZEnum || name==MaterialsRheologyZbarEnum){
    1694                                 matice->inputs->AddInput(new TriaP1Input(name,values));
     1694                                material->inputs->AddInput(new TriaP1Input(name,values));
    16951695                        }
    16961696                        else{
     
    27922792        switch(response_enum){
    27932793                case MaterialsRheologyBbarEnum:
    2794                         *presponse=this->matice->GetBbar();
     2794                        *presponse=this->material->GetBbar();
    27952795                        break;
    27962796                case MaterialsRheologyZbarEnum:
    2797                         *presponse=this->matice->GetZbar();
     2797                        *presponse=this->material->GetZbar();
    27982798                        break;
    27992799                case VelEnum:{
     
    29092909                this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    29102910                this->GetStrainRate2d(&oldepsilon[0],&xyz_list[0][0],gauss,vxold_input,vyold_input);
    2911                 matice->GetViscosity2d(&viscosity, &epsilon[0]);
    2912                 matice->GetViscosity2d(&oldviscosity, &oldepsilon[0]);
     2911                material->GetViscosity2d(&viscosity, &epsilon[0]);
     2912                material->GetViscosity2d(&oldviscosity, &oldepsilon[0]);
    29132913                thickness_input->GetInputValue(&thickness, gauss);
    29142914
     
    30923092        rho_ice=matpar->GetRhoIce();
    30933093        gravity=matpar->GetG();
    3094         n=matice->GetN();
    3095         B=matice->GetBbar();
     3094        n=material->GetN();
     3095        B=material->GetBbar();
    30963096        Input* slopex_input=inputs->GetInput(SurfaceSlopeXEnum); _assert_(slopex_input);
    30973097        Input* slopey_input=inputs->GetInput(SurfaceSlopeYEnum); _assert_(slopey_input);
     
    31643164                thickness_input->GetInputValue(&thickness, gauss);
    31653165                this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    3166                 matice->GetViscosity2dDerivativeEpsSquare(&mu_prime,&epsilon[0]);
     3166                material->GetViscosity2dDerivativeEpsSquare(&mu_prime,&epsilon[0]);
    31673167                eps1[0]=2*epsilon[0]+epsilon[1];   eps2[0]=epsilon[2];
    31683168                eps1[1]=epsilon[2];                eps2[1]=epsilon[0]+2*epsilon[1];
     
    34093409
    34103410                if(control_type[i]==MaterialsRheologyBbarEnum || control_type[i]==MaterialsRheologyZbarEnum){
    3411                         input=(Input*)matice->inputs->GetInput(control_type[i]); _assert_(input);
     3411                        input=(Input*)material->inputs->GetInput(control_type[i]); _assert_(input);
    34123412                }
    34133413                else{
     
    34363436
    34373437        if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum){
    3438                 input=(Input*)matice->inputs->GetInput(enum_type);
     3438                input=(Input*)material->inputs->GetInput(enum_type);
    34393439        }
    34403440        else{
     
    34543454
    34553455        if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum){
    3456                 input=(Input*)matice->inputs->GetInput(enum_type);
     3456                input=(Input*)material->inputs->GetInput(enum_type);
    34573457        }
    34583458        else{
     
    34733473
    34743474        if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum){
    3475                 input=(Input*)matice->inputs->GetInput(enum_type);
     3475                input=(Input*)material->inputs->GetInput(enum_type);
    34763476        }
    34773477        else{
     
    35673567        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    35683568        GradientIndexing(&doflist1[0],control_index);
    3569         Input* rheologyb_input=matice->inputs->GetInput(MaterialsRheologyBbarEnum); _assert_(rheologyb_input);
     3569        Input* rheologyb_input=material->inputs->GetInput(MaterialsRheologyBbarEnum); _assert_(rheologyb_input);
    35703570        Input* weights_input=inputs->GetInput(InversionCostFunctionsCoefficientsEnum);                _assert_(weights_input);
    35713571
     
    36073607        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    36083608        GradientIndexing(&doflist1[0],control_index);
    3609         Input* rheologyz_input=matice->inputs->GetInput(MaterialsRheologyZbarEnum); _assert_(rheologyz_input);
     3609        Input* rheologyz_input=material->inputs->GetInput(MaterialsRheologyZbarEnum); _assert_(rheologyz_input);
    36103610        Input* weights_input=inputs->GetInput(InversionCostFunctionsCoefficientsEnum);                _assert_(weights_input);
    36113611
     
    36563656        Input* adjointx_input=inputs->GetInput(AdjointxEnum);                       _assert_(adjointx_input);
    36573657        Input* adjointy_input=inputs->GetInput(AdjointyEnum);                       _assert_(adjointy_input);
    3658         Input* rheologyb_input=matice->inputs->GetInput(MaterialsRheologyBbarEnum); _assert_(rheologyb_input);
     3658        Input* rheologyb_input=material->inputs->GetInput(MaterialsRheologyBbarEnum); _assert_(rheologyb_input);
    36593659
    36603660        /* Start  looping on the number of gaussian points: */
     
    36723672
    36733673                this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    3674                 matice->GetViscosityComplement(&viscosity_complement,&epsilon[0]);
     3674                material->GetViscosityComplement(&viscosity_complement,&epsilon[0]);
    36753675
    36763676                GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0],gauss);
     
    37133713        Input* adjointx_input=inputs->GetInput(AdjointxEnum);                       _assert_(adjointx_input);
    37143714        Input* adjointy_input=inputs->GetInput(AdjointyEnum);                       _assert_(adjointy_input);
    3715         Input* rheologyz_input=matice->inputs->GetInput(MaterialsRheologyZbarEnum); _assert_(rheologyz_input);
     3715        Input* rheologyz_input=material->inputs->GetInput(MaterialsRheologyZbarEnum); _assert_(rheologyz_input);
    37163716
    37173717        /* Start  looping on the number of gaussian points: */
     
    37293729
    37303730                this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    3731                 matice->GetViscosityZComplement(&viscosity_complement,&epsilon[0]);
     3731                material->GetViscosityZComplement(&viscosity_complement,&epsilon[0]);
    37323732
    37333733                GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0],gauss);
     
    40134013        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    40144014        Input* weights_input  =inputs->GetInput(InversionCostFunctionsCoefficientsEnum);              _assert_(weights_input);
    4015         Input* rheologyb_input=matice->inputs->GetInput(MaterialsRheologyBbarEnum); _assert_(rheologyb_input);
     4015        Input* rheologyb_input=material->inputs->GetInput(MaterialsRheologyBbarEnum); _assert_(rheologyb_input);
    40164016
    40174017        /* Start looping on the number of gaussian points: */
     
    51105110                thickness_input->GetInputValue(&thickness, gauss);
    51115111                this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    5112                 matice->GetViscosity2dDerivativeEpsSquare(&mu_prime,&epsilon[0]);
     5112                material->GetViscosity2dDerivativeEpsSquare(&mu_prime,&epsilon[0]);
    51135113                eps1[0]=2*epsilon[0]+epsilon[1];   eps2[0]=epsilon[2];
    51145114                eps1[1]=epsilon[2];                eps2[1]=epsilon[0]+2*epsilon[1];
     
    52165216        /*Get input (either in element or material)*/
    52175217        if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheologyZbarEnum){
    5218                 input=(Input*)matice->inputs->GetInput(control_enum); _assert_(input);
     5218                input=(Input*)material->inputs->GetInput(control_enum); _assert_(input);
    52195219        }
    52205220        else{
     
    52515251
    52525252        if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheologyZbarEnum){
    5253                 input=(Input*)matice->inputs->GetInput(control_enum); _assert_(input);
     5253                input=(Input*)material->inputs->GetInput(control_enum); _assert_(input);
    52545254        }
    52555255        else{
  • issm/trunk-jpl/src/c/classes/objects/Elements/Tria.h

    r13119 r13129  
    1515class IoModel;
    1616class Node;
    17 class Matice;
     17class Material;
    1818class Matpar;
    1919class ElementMatrix;
     
    3333
    3434                Node   **nodes;    // 3 nodes
    35                 Matice  *matice;   // 1 material ice
     35                Material  *material;   // 1 material ice
    3636                Matpar  *matpar;   // 1 material parameter
    3737                int      horizontalneighborsids[3];
  • issm/trunk-jpl/src/c/classes/objects/Elements/TriaHook.cpp

    r12832 r13129  
    2525        numanalyses=UNDEF;
    2626        this->hnodes=NULL;
    27         this->hmatice=NULL;
     27        this->hmaterial=NULL;
    2828        this->hmatpar=NULL;
    2929}
     
    3737        }
    3838        delete [] this->hnodes;
    39         delete hmatice;
     39        delete hmaterial;
    4040        delete hmatpar;
    4141
    4242}
    4343/*}}}*/
    44 /*FUNCTION TriaHook::TriaHook(int in_numanalyses,int matice_id, int matpar_id){{{*/
    45 TriaHook::TriaHook(int in_numanalyses,int matice_id, IoModel* iomodel){
     44/*FUNCTION TriaHook::TriaHook(int in_numanalyses,int material_id, int matpar_id){{{*/
     45TriaHook::TriaHook(int in_numanalyses,int material_id, IoModel* iomodel){
    4646
    4747        /*intermediary: */
     
    5353        this->numanalyses=in_numanalyses;
    5454        this->hnodes= new Hook*[in_numanalyses];
    55         this->hmatice=new Hook(&matice_id,1);
     55        this->hmaterial=new Hook(&material_id,1);
    5656        this->hmatpar=new Hook(&matpar_id,1);
    5757
  • issm/trunk-jpl/src/c/classes/objects/Elements/TriaHook.h

    r12365 r13129  
    1414                int    numanalyses; //number of analysis types
    1515                Hook** hnodes; // 3 nodes for each analysis type
    16                 Hook*  hmatice; // 1 ice material
     16                Hook*  hmaterial; // 1 ice material
    1717                Hook*  hmatpar; // 1 material parameter
    1818
     
    2020                /*FUNCTION constructors, destructors {{{*/
    2121                TriaHook();
    22                 TriaHook(int in_numanalyses,int matice_id, IoModel* iomodel);
     22                TriaHook(int in_numanalyses,int material_id, IoModel* iomodel);
    2323                ~TriaHook();
    2424                void SetHookNodes(int* node_ids,int analysis_counter);
  • issm/trunk-jpl/src/c/classes/objects/Loads/Icefront.cpp

    r13036 r13129  
    498498
    499499        /*clean-up and return*/
    500         delete tria->matice;
     500        delete tria->material;
    501501        delete tria;
    502502        delete icefront;
  • issm/trunk-jpl/src/c/classes/objects/Materials/Material.h

    r12832 r13129  
    1717
    1818        public:
     19                Inputs*  inputs;
    1920                virtual       ~Material(){};
     21                /*WARNING: input should not be public but it is an easy way to update B from T (using UpdateFromSolution) from Pentas*/
    2022
    2123                /*Numerics*/
    22                 virtual void   InputDuplicate(int original_enum,int new_enum)=0;
    23                 virtual void   Configure(Elements* elements)=0;
    24                 virtual void   GetVectorFromInputs(Vector* vector,int input_enum)=0;
     24                virtual void       InputDuplicate(int original_enum,int new_enum)=0;
     25                virtual void       Configure(Elements* elements)=0;
     26                virtual void       GetVectorFromInputs(Vector* vector,int input_enum)=0;
     27                virtual void       GetViscosity2d(IssmDouble* pviscosity, IssmDouble* pepsilon)=0;
     28                virtual void       GetViscosity3d(IssmDouble* pviscosity3d, IssmDouble* pepsilon)=0;
     29                virtual void       GetViscosity3dStokes(IssmDouble* pviscosity3d, IssmDouble* epsilon)=0;
     30                virtual void       GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon)=0;
     31                virtual void       GetViscosityZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon)=0;
     32                virtual void       GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon)=0;
     33                virtual void       GetViscosity2dDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon)=0;
     34                virtual IssmDouble GetA()=0;
     35                virtual IssmDouble GetB()=0;
     36                virtual IssmDouble GetBbar()=0;
     37                virtual IssmDouble GetN()=0;
     38                virtual IssmDouble GetZ()=0;
     39                virtual IssmDouble GetZbar()=0;
    2540
    2641};
  • issm/trunk-jpl/src/c/classes/objects/Materials/Matice.cpp

    r13119 r13129  
    172172}
    173173/*}}}*/
    174 /*FUNCTION Matice::GetZ {{{*/
    175 IssmDouble Matice::GetZ(){
    176 
    177         /*Output*/
    178         IssmDouble Z;
    179 
    180         inputs->GetInputAverage(&Z,MaterialsRheologyZEnum);
    181         return Z;
    182 }
    183 /*}}}*/
    184 /*FUNCTION Matice::GetZbar {{{*/
    185 IssmDouble Matice::GetZbar(){
    186 
    187         /*Output*/
    188         IssmDouble Zbar;
    189 
    190         inputs->GetInputAverage(&Zbar,MaterialsRheologyZbarEnum);
    191         return Zbar;
    192 }
    193 /*}}}*/
    194174/*FUNCTION Matice::GetVectorFromInputs{{{*/
    195175void  Matice::GetVectorFromInputs(Vector* vector,int input_enum){
     
    227207void  Matice::GetViscosity2d(IssmDouble* pviscosity, IssmDouble* epsilon){
    228208        /*From a string tensor and a material object, return viscosity, using Glen's flow law.
    229                                                                                                    Z * B
     209                                                                                                    B
    230210          viscosity= -------------------------------------------------------------------
    231211                                                  2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
     
    246226        /*Intermediary: */
    247227        IssmDouble A,e;
    248         IssmDouble Btmp,B,n,Z;
     228        IssmDouble B,n;
    249229
    250230        /*Get B and n*/
    251         Btmp=GetBbar();
    252         Z=GetZbar();
     231        B=GetBbar();
    253232        n=GetN();
    254         B=Z*Btmp;
    255233
    256234        if (n==1){
     
    314292        /*Intermediaries: */
    315293        IssmDouble A,e;
    316         IssmDouble B,n,Z;
    317 
    318         /*Get B, Z and n*/
     294        IssmDouble B,n;
     295
     296        /*Get B and n*/
     297        B=GetB();
    319298        n=GetN();
    320         Z=GetZ();
    321         B=Z*GetB();
    322299
    323300        if (n==1){
     
    385362        /*Intermediaries: */
    386363        IssmDouble A,e;
    387         IssmDouble B,n,Z;
     364        IssmDouble B,n;
    388365        IssmDouble eps0;
    389366
    390367        /*Get B and n*/
    391368        eps0=pow((IssmDouble)10,(IssmDouble)-27);
     369        B=GetB();
    392370        n=GetN();
    393         Z=GetZ();
    394         B=Z*GetB();
    395371       
    396372        if (n==1){
     
    476452               
    477453                        viscosity_complement=1/(2*pow(A,e));
    478                 }
    479         }
    480         else{
    481                 viscosity_complement=4.5*pow((IssmDouble)10,(IssmDouble)17);
    482         }
    483 
    484         /*Checks in debugging mode*/
    485         _assert_(B>0);
    486         _assert_(n>0);
    487         _assert_(viscosity_complement>0);
    488                
    489         /*Return: */
    490         *pviscosity_complement=viscosity_complement;
    491 }
    492 /*}}}*/
    493 /*FUNCTION Matice::GetViscosityZComplement {{{*/
    494 void  Matice::GetViscosityZComplement(IssmDouble* pviscosity_complement, IssmDouble* epsilon){
    495         /*Return viscosity derivative for control method d(mu)/dZ:
    496          *
    497          *                                                                                             B
    498          * dviscosity= -------------------------------------------------------------------
    499          *                                2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
    500          *
    501          * If epsilon is NULL, it means this is the first time Gradjb is being run, and we
    502          * return mu20, initial viscosity.
    503          */
    504        
    505         /*output: */
    506         IssmDouble viscosity_complement;
    507 
    508         /*input strain rate: */
    509         IssmDouble exx,eyy,exy;
    510 
    511         /*Intermediary value A and exponent e: */
    512         IssmDouble A,e;
    513         IssmDouble B,n;
    514 
    515         /*Get B and n*/
    516         B=GetBbar();
    517         n=GetN();
    518 
    519         if(epsilon){
    520                 exx=*(epsilon+0);
    521                 eyy=*(epsilon+1);
    522                 exy=*(epsilon+2);
    523 
    524                 /*Build viscosity: mu2=B/(2*A^e) */
    525                 A=pow(exx,2)+pow(eyy,2)+pow(exy,2)+exx*eyy;
    526                 if(A==0){
    527                         /*Maximum viscosity_complement for 0 shear areas: */
    528                         viscosity_complement=2.25*pow((IssmDouble)10,(IssmDouble)17);
    529                 }
    530                 else{
    531                         e=(n-1)/(2*n);
    532                
    533                         viscosity_complement=B/(2*pow(A,e));
    534454                }
    535455        }
     
    781701                }
    782702
    783                 /*Get Z*/
    784                 if (iomodel->Data(MaterialsRheologyZEnum)) {
    785                         for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheologyZEnum)[reCast<int,IssmDouble>(iomodel->Data(MeshElementsEnum)[num_vertices*index+i]-1)];
    786                         this->inputs->AddInput(new TriaP1Input(MaterialsRheologyZbarEnum,nodeinputs));
    787                 }
    788 
    789703                /*Control Inputs*/
    790704                #ifdef _HAVE_CONTROL_
     
    801715                                                }
    802716                                                break;
    803                                         case MaterialsRheologyZbarEnum:
    804                                                 if (iomodel->Data(MaterialsRheologyZEnum)){
    805                                                         _assert_(iomodel->Data(MaterialsRheologyZEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));
    806                                                         for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(MaterialsRheologyZEnum)[int(iomodel->Data(MeshElementsEnum)[num_vertices*index+j]-1)];
    807                                                         for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->Data(InversionMinParametersEnum)[int(iomodel->Data(MeshElementsEnum)[num_vertices*index+j]-1)*num_control_type+i];
    808                                                         for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->Data(InversionMaxParametersEnum)[int(iomodel->Data(MeshElementsEnum)[num_vertices*index+j]-1)*num_control_type+i];
    809                                                         this->inputs->AddInput(new ControlInput(MaterialsRheologyZbarEnum,TriaP1InputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    810                                                 }
    811                                                 break;
    812 
    813717                                }
    814718                        }
     
    837741                        for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheologyNEnum)[index];
    838742                        this->inputs->AddInput(new PentaP1Input(MaterialsRheologyNEnum,nodeinputs));
    839                 }
    840 
    841                 /*Get Z*/
    842                 if (iomodel->Data(MaterialsRheologyZEnum)) {
    843                         for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheologyZEnum)[int(iomodel->Data(MeshElementsEnum)[num_vertices*index+i]-1)];
    844                         this->inputs->AddInput(new PentaP1Input(MaterialsRheologyZEnum,nodeinputs));
    845743                }
    846744
     
    859757                                                }
    860758                                                break;
    861                                         case MaterialsRheologyZbarEnum:
    862                                                 if (iomodel->Data(MaterialsRheologyZEnum)){
    863                                                         _assert_(iomodel->Data(MaterialsRheologyZEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));
    864                                                         for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(MaterialsRheologyZEnum)[int(iomodel->Data(MeshElementsEnum)[num_vertices*index+j]-1)];
    865                                                         for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->Data(InversionMinParametersEnum)[int(iomodel->Data(MeshElementsEnum)[num_vertices*index+j]-1)*num_control_type+i];
    866                                                         for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->Data(InversionMaxParametersEnum)[int(iomodel->Data(MeshElementsEnum)[num_vertices*index+j]-1)*num_control_type+i];
    867                                                         this->inputs->AddInput(new ControlInput(MaterialsRheologyZEnum,PentaP1InputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    868                                                 }
    869                                                 break;
    870759                                }
    871760                        }
     
    886775                                name==MaterialsRheologyBEnum ||
    887776                                name==MaterialsRheologyBbarEnum ||
    888                                 name==MaterialsRheologyNEnum ||
    889                                 name==MaterialsRheologyZEnum ||
    890                                 name==MaterialsRheologyZbarEnum
     777                                name==MaterialsRheologyNEnum
    891778                ){
    892779                return true;
  • issm/trunk-jpl/src/c/classes/objects/Materials/Matice.h

    r13119 r13129  
    1515
    1616        private:
    17                 /*Id*/
    1817                int        mid;
    19 
    20                 /*hooks: */
    2118                Hook* helement;
    2219
    2320        public:
    24                 /*WARNING: input should not be public but it is an easy way to update B from T (using UpdateFromSolution) from Pentas*/
    25                 Inputs*  inputs;
    26 
    2721                /*Matice constructors, destructors: {{{*/
    2822                Matice();
     
    5650                void   Configure(Elements* elements);
    5751                void   GetVectorFromInputs(Vector* vector,int input_enum);
    58                 /*}}}*/
    59                 /*Matice Numerics: {{{*/
    60                 void   SetCurrentConfiguration(Elements* elementsin,Loads* loadsin,Nodes* nodesin,Vertices* verticesin,Materials* materialsin,Parameters* parametersin);
    61                 void   GetViscosity2d(IssmDouble* pviscosity, IssmDouble* pepsilon);
    62                 void   GetViscosity3d(IssmDouble* pviscosity3d, IssmDouble* pepsilon);
    63                 void   GetViscosity3dStokes(IssmDouble* pviscosity3d, IssmDouble* epsilon);
    64                 void   GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
    65                 void   GetViscosityZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
    66                 void   GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
    67                 void   GetViscosity2dDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
     52                void       SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
     53                void       GetViscosity2d(IssmDouble* pviscosity, IssmDouble* pepsilon);
     54                void       GetViscosity3d(IssmDouble* pviscosity3d, IssmDouble* pepsilon);
     55                void       GetViscosity3dStokes(IssmDouble* pviscosity3d, IssmDouble* epsilon);
     56                void       GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
     57                void GetViscosityZComplement(IssmDouble*, IssmDouble*){_error_("not supported");};
     58                void       GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
     59                void       GetViscosity2dDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
    6860                IssmDouble GetA();
    6961                IssmDouble GetB();
    7062                IssmDouble GetBbar();
     63                IssmDouble GetZ(){_error_("not supported");};
     64                IssmDouble GetZbar(){_error_("not supported");};
    7165                IssmDouble GetN();
    72                 IssmDouble GetZ();
    73                 IssmDouble GetZbar();
    74                 bool   IsInput(int name);
     66                bool       IsInput(int name);
    7567                /*}}}*/
    7668};
  • issm/trunk-jpl/src/c/classes/objects/Materials/Matpar.cpp

    r12832 r13129  
    4545        iomodel->Constant(&this->hydro_p,HydrologyPEnum);
    4646        iomodel->Constant(&this->hydro_q,HydrologyQEnum);
     47        this->inputs=NULL; /*not used here*/
    4748}
    4849/*}}}*/
  • issm/trunk-jpl/src/c/classes/objects/Materials/Matpar.h

    r13056 r13129  
    6868                void   Configure(Elements* elements);
    6969                void   GetVectorFromInputs(Vector* vector,int input_enum){return;}
     70                void       GetViscosity2d(IssmDouble* pviscosity, IssmDouble* pepsilon){_error_("not supported");};
     71                void       GetViscosity3d(IssmDouble* pviscosity3d, IssmDouble* pepsilon){_error_("not supported");};
     72                void       GetViscosity3dStokes(IssmDouble* pviscosity3d, IssmDouble* epsilon){_error_("not supported");};
     73                void       GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon){_error_("not supported");};
     74                void       GetViscosityZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon){_error_("not supported");};
     75                void       GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon){_error_("not supported");};
     76                void       GetViscosity2dDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon){_error_("not supported");};
     77                IssmDouble GetA(){_error_("not supported");};
     78                IssmDouble GetB(){_error_("not supported");};
     79                IssmDouble GetBbar(){_error_("not supported");};
     80                IssmDouble GetN(){_error_("not supported");};
     81                IssmDouble GetZ(){_error_("not supported");};
     82                IssmDouble GetZbar(){_error_("not supported");};
    7083                /*}}}*/
    7184                /*Numerics: {{{*/
  • issm/trunk-jpl/src/c/classes/objects/Params/BoolParam.h

    r13036 r13129  
    4343                int   InstanceEnum(){return enum_type;}
    4444                void  GetParameterValue(bool* pbool){*pbool=value;}
    45                 void  GetParameterValue(int* pinteger){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return an integer");}
    46                 void  GetParameterValue(int** pintarray,int* pM){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return an array of integers");}
    47                 void  GetParameterValue(int** pintarray,int* pM,int* pN){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return an array of integers");}
    48                 void  GetParameterValue(IssmDouble* pIssmDouble){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a IssmDouble");}
     45                void  GetParameterValue(int* pinteger){_error_("Param "<< EnumToStringx(enum_type) << " cannot return an integer");}
     46                void  GetParameterValue(int** pintarray,int* pM){_error_("Param "<< EnumToStringx(enum_type) << " cannot return an array of integers");}
     47                void  GetParameterValue(int** pintarray,int* pM,int* pN){_error_("Param "<< EnumToStringx(enum_type) << " cannot return an array of integers");}
     48                void  GetParameterValue(IssmDouble* pIssmDouble){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a IssmDouble");}
    4949                void  GetParameterValue(IssmDouble* pdouble,IssmDouble time){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a IssmDouble for a given time");}
    50                 void  GetParameterValue(char** pstring){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a string");}
    51                 void  GetParameterValue(char*** pstringarray,int* pM){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a string array");}
    52                 void  GetParameterValue(IssmDouble** pIssmDoublearray,int* pM){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a IssmDouble array");}
    53                 void  GetParameterValue(IssmDouble** pIssmDoublearray,int* pM, int* pN){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a IssmDouble array");}
    54                 void  GetParameterValue(IssmDouble*** parray, int* pM,int** pmdims, int** pndims){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a matrix array");}
    55                 void  GetParameterValue(Vector** pvec){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a Vec");}
    56                 void  GetParameterValue(Matrix** pmat){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a Mat");}
    57                 void  GetParameterValue(FILE** pfid){_error_("Param "<< EnumToStringx(enum_type) << ") cannot return a FILE");}
     50                void  GetParameterValue(char** pstring){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a string");}
     51                void  GetParameterValue(char*** pstringarray,int* pM){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a string array");}
     52                void  GetParameterValue(IssmDouble** pIssmDoublearray,int* pM){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a IssmDouble array");}
     53                void  GetParameterValue(IssmDouble** pIssmDoublearray,int* pM, int* pN){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a IssmDouble array");}
     54                void  GetParameterValue(IssmDouble*** parray, int* pM,int** pmdims, int** pndims){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a matrix array");}
     55                void  GetParameterValue(Vector** pvec){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a Vec");}
     56                void  GetParameterValue(Matrix** pmat){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a Mat");}
     57                void  GetParameterValue(FILE** pfid){_error_("Param "<< EnumToStringx(enum_type) << " cannot return a FILE");}
    5858
    5959                void  SetValue(bool boolean){this->value=boolean;}
    60                 void  SetValue(int integer){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold an int");}
    61                 void  SetValue(IssmDouble scalar){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold an IssmPDouble");}
    62                 void  SetValue(char* string){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a string");}
    63                 void  SetValue(char** stringarray,int M){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a string array");}
    64                 void  SetValue(IssmDouble* IssmDoublearray,int M){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a IssmDouble array");}
    65                 void  SetValue(IssmDouble* pIssmDoublearray,int M,int N){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a IssmDouble array");}
    66                 void  SetValue(int* intarray,int M){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a int array");}
    67                 void  SetValue(int* pintarray,int M,int N){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a int array");}
    68                 void  SetValue(Vector* vec){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a Vec");}
    69                 void  SetValue(Matrix* mat){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a Mat");}
    70                 void  SetValue(FILE* fid){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold a FILE");}
    71                 void  SetValue(IssmDouble** array, int M, int* mdim_array, int* ndim_array){_error_("Param "<< EnumToStringx(enum_type) << ") cannot hold an array of matrices");}
     60                void  SetValue(int integer){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold an int");}
     61                void  SetValue(IssmDouble scalar){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold an IssmPDouble");}
     62                void  SetValue(char* string){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a string");}
     63                void  SetValue(char** stringarray,int M){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a string array");}
     64                void  SetValue(IssmDouble* IssmDoublearray,int M){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a IssmDouble array");}
     65                void  SetValue(IssmDouble* pIssmDoublearray,int M,int N){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a IssmDouble array");}
     66                void  SetValue(int* intarray,int M){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a int array");}
     67                void  SetValue(int* pintarray,int M,int N){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a int array");}
     68                void  SetValue(Vector* vec){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a Vec");}
     69                void  SetValue(Matrix* mat){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a Mat");}
     70                void  SetValue(FILE* fid){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold a FILE");}
     71                void  SetValue(IssmDouble** array, int M, int* mdim_array, int* ndim_array){_error_("Param "<< EnumToStringx(enum_type) << " cannot hold an array of matrices");}
    7272                void  UnitConversion(int direction_enum);
    7373               
  • issm/trunk-jpl/src/c/classes/objects/objects.h

    r12835 r13129  
    108108#include "./Materials/Material.h"
    109109#include "./Materials/Matice.h"
     110#include "./Materials/Matdamageice.h"
    110111#include "./Materials/Matpar.h"
    111112
  • issm/trunk-jpl/src/c/modules/EnumToStringx/EnumToStringx.cpp

    r13119 r13129  
    279279                case MacAyeal3dIceFrontEnum : return "MacAyeal3dIceFront";
    280280                case MaticeEnum : return "Matice";
     281                case MatdamageiceEnum : return "Matdamageice";
    281282                case MatparEnum : return "Matpar";
    282283                case NodeEnum : return "Node";
  • issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateElementsVerticesAndMaterials.cpp

    r13119 r13129  
    1717        /*Intermediary*/
    1818        int i,j,k,n;
    19         int dim;
     19        int dim,materials_type;
    2020        int numberofelements;
    2121        int numberofvertices;
     
    2323
    2424        /*DataSets: */
    25         Elements*     elements  = NULL;
    26         Vertices*     vertices = NULL;
    27         Materials*    materials = NULL;
     25        Elements  *elements  = NULL;
     26        Vertices  *vertices = NULL;
     27        Materials *materials = NULL;
    2828
    2929        /*Fetch parameters: */
     
    3232        iomodel->Constant(&numberofvertices,MeshNumberofverticesEnum);
    3333        iomodel->Constant(&control_analysis,InversionIscontrolEnum);
     34        iomodel->Constant(&materials_type,MaterialsEnum);
    3435
    3536        /*Did we already create the elements? : */
     
    4445        ElementsAndVerticesPartitioning(&iomodel->my_elements,&iomodel->my_vertices,iomodel);
    4546       
    46         /*Fetch data needed: */
    47         iomodel->FetchData(5,MeshElementsEnum,MeshElementconnectivityEnum,MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheologyZEnum);
     47        iomodel->FetchData(2,MeshElementsEnum,MeshElementconnectivityEnum);
    4848        #ifdef _HAVE_3D_
    4949        if(dim==3)iomodel->FetchData(2,MeshUpperelementsEnum,MeshLowerelementsEnum);
     
    5151        if(control_analysis)iomodel->FetchData(3,InversionControlParametersEnum,InversionMinParametersEnum,InversionMaxParametersEnum);
    5252       
    53         /*Create elements and materials: */
     53        /*Create elements*/
    5454        for (i=0;i<numberofelements;i++){
    5555                if(iomodel->my_elements[i]){
     
    6060                        else       elements->AddObject(new Penta(i+1,i,i,iomodel,nummodels));
    6161                #endif
    62 
    63                         /*Create and add material property to materials dataset: */
    64                         materials->AddObject(new Matice(i+1,i,iomodel));
    6562                }
    6663        }
    6764       
     65        /*Create materials*/
     66        switch(materials_type){
     67                case MaticeEnum:
     68                        iomodel->FetchData(2,MaterialsRheologyBEnum,MaterialsRheologyNEnum);
     69                        for (i=0;i<numberofelements;i++) if(iomodel->my_elements[i]) materials->AddObject(new Matice(i+1,i,iomodel));
     70                        break;
     71                case MatdamageiceEnum:
     72                        iomodel->FetchData(3,MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheologyZEnum);
     73                        for (i=0;i<numberofelements;i++) if(iomodel->my_elements[i]) materials->AddObject(new Matdamageice(i+1,i,iomodel));
     74                        break;
     75                default:
     76                        _error_("Materials "<<EnumToStringx(materials_type)<<" not supported");
     77        }
     78
    6879        /*Free data: */
    6980        iomodel->DeleteData(10,MeshElementsEnum,MeshElementconnectivityEnum,MeshUpperelementsEnum,MeshLowerelementsEnum,
  • issm/trunk-jpl/src/c/modules/ModelProcessorx/DiagnosticHoriz/UpdateElementsDiagnosticHoriz.cpp

    r13119 r13129  
    1616void    UpdateElementsDiagnosticHoriz(Elements* elements, IoModel* iomodel,int analysis_counter,int analysis_type){
    1717
    18         int    dim;
     18        int    dim,materials_type;
    1919        int    numberofelements;
    2020        bool   ismacayealpattyn;
     
    3232        iomodel->Constant(&control_analysis,InversionIscontrolEnum);
    3333        iomodel->Constant(&dakota_analysis,QmuIsdakotaEnum);
     34        iomodel->Constant(&materials_type,MaterialsEnum);
    3435
    3536        /*Now, is the flag macayaealpattyn on? otherwise, do nothing: */
     
    6162        iomodel->FetchDataToInput(elements,MaterialsRheologyBEnum);
    6263        iomodel->FetchDataToInput(elements,MaterialsRheologyNEnum);
    63         iomodel->FetchDataToInput(elements,MaterialsRheologyZEnum);
    6464        iomodel->FetchDataToInput(elements,VxEnum);
    6565        iomodel->FetchDataToInput(elements,VyEnum);
    66 
     66        if(materials_type==MatdamageiceEnum){
     67                iomodel->FetchDataToInput(elements,MaterialsRheologyZEnum);
     68        }
    6769        if (dim==3){
    6870                iomodel->FetchDataToInput(elements,MeshElementonbedEnum);
  • issm/trunk-jpl/src/c/modules/StringToEnumx/StringToEnumx.cpp

    r13119 r13129  
    286286              else if (strcmp(name,"MacAyeal3dIceFront")==0) return MacAyeal3dIceFrontEnum;
    287287              else if (strcmp(name,"Matice")==0) return MaticeEnum;
     288              else if (strcmp(name,"Matdamageice")==0) return MatdamageiceEnum;
    288289              else if (strcmp(name,"Matpar")==0) return MatparEnum;
    289290              else if (strcmp(name,"Node")==0) return NodeEnum;
     
    383384              else if (strcmp(name,"Enthalpy")==0) return EnthalpyEnum;
    384385              else if (strcmp(name,"EnthalpyPicard")==0) return EnthalpyPicardEnum;
    385               else if (strcmp(name,"ThicknessAbsGradient")==0) return ThicknessAbsGradientEnum;
    386386         else stage=4;
    387387   }
    388388   if(stage==4){
    389               if (strcmp(name,"ThicknessAlongGradient")==0) return ThicknessAlongGradientEnum;
     389              if (strcmp(name,"ThicknessAbsGradient")==0) return ThicknessAbsGradientEnum;
     390              else if (strcmp(name,"ThicknessAlongGradient")==0) return ThicknessAlongGradientEnum;
    390391              else if (strcmp(name,"ThicknessAcrossGradient")==0) return ThicknessAcrossGradientEnum;
    391392              else if (strcmp(name,"StepResponses")==0) return StepResponsesEnum;
Note: See TracChangeset for help on using the changeset viewer.