Changeset 16146


Ignore:
Timestamp:
09/17/13 12:44:19 (12 years ago)
Author:
Eric.Larour
Message:

CHG: switched damage representation from Z to D=(1-Z). Z is now 1/(1-D) instead of the original (1-D).
which was incorrectly (on purpose) named.

Location:
issm/trunk-jpl
Files:
2 added
25 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

    r16145 r16146  
    10061006        /*Computeportion of the element that is grounded*/
    10071007
    1008         int         normal_orientation=0;
     1008        int         normal_orientation;
    10091009        IssmDouble  s1,s2;
    10101010        IssmDouble  levelset[NUMVERTICES];
     
    10131013        GetInputListOnVertices(&levelset[0],levelsetenum);
    10141014
    1015         if(levelset[0]*levelset[1]>0.){ //Nodes 0 and 1 are similar, so points must be found on segment 0-2 and 1-2
     1015        if(levelset[0]*levelset[1]>0){ //Nodes 0 and 1 are similar, so points must be found on segment 0-2 and 1-2
    10161016                /*Portion of the segments*/
    10171017                s1=levelset[2]/(levelset[2]-levelset[1]);
    10181018                s2=levelset[2]/(levelset[2]-levelset[0]);
    10191019
    1020                 if(levelset[2]>0.) normal_orientation=1;
     1020                if(levelset[2]>0) normal_orientation=1;
    10211021                /*New point 1*/
    10221022                xyz_zero[3*normal_orientation+0]=xyz_list[2][0]+s1*(xyz_list[1][0]-xyz_list[2][0]);
     
    10291029                xyz_zero[3*(1-normal_orientation)+2]=xyz_list[2][2]+s2*(xyz_list[0][2]-xyz_list[2][2]);
    10301030        }
    1031         else if(levelset[1]*levelset[2]>0.){ //Nodes 1 and 2 are similar, so points must be found on segment 0-1 and 0-2
     1031        else if(levelset[1]*levelset[2]>0){ //Nodes 1 and 2 are similar, so points must be found on segment 0-1 and 0-2
    10321032                /*Portion of the segments*/
    10331033                s1=levelset[0]/(levelset[0]-levelset[2]);
    10341034                s2=levelset[0]/(levelset[0]-levelset[1]);
    10351035
    1036                 if(levelset[0]>0.) normal_orientation=1;
     1036                if(levelset[0]>0) normal_orientation=1;
    10371037                /*New point 1*/
    10381038                xyz_zero[3*normal_orientation+0]=xyz_list[0][0]+s1*(xyz_list[2][0]-xyz_list[0][0]);
     
    10451045                xyz_zero[3*(1-normal_orientation)+2]=xyz_list[0][2]+s2*(xyz_list[1][2]-xyz_list[0][2]);
    10461046        }
    1047         else if(levelset[0]*levelset[2]>0.){ //Nodes 0 and 2 are similar, so points must be found on segment 1-0 and 1-2
     1047        else if(levelset[0]*levelset[2]>0){ //Nodes 0 and 2 are similar, so points must be found on segment 1-0 and 1-2
    10481048                /*Portion of the segments*/
    10491049                s1=levelset[1]/(levelset[1]-levelset[0]);
    10501050                s2=levelset[1]/(levelset[1]-levelset[2]);
    10511051
    1052                 if(levelset[1]>0.) normal_orientation=1;
     1052                if(levelset[1]>0) normal_orientation=1;
    10531053                /*New point 0*/
    10541054                xyz_zero[3*normal_orientation+0]=xyz_list[1][0]+s1*(xyz_list[0][0]-xyz_list[1][0]);
     
    10611061                xyz_zero[3*(1-normal_orientation)+2]=xyz_list[1][2]+s2*(xyz_list[2][2]-xyz_list[1][2]);
    10621062        }
    1063         else if(levelset[0]==0. && levelset[1]==0.){ //front is on point 0 and 1
     1063        else if(levelset[0]==0 && levelset[1]==0){ //front is on point 0 and 1
    10641064                xyz_zero[3*0+0]=xyz_list[0][0];
    10651065                xyz_zero[3*0+1]=xyz_list[0][1];
     
    10711071                xyz_zero[3*1+2]=xyz_list[1][2];
    10721072        }
    1073         else if(levelset[0]==0. && levelset[2]==0.){ //front is on point 0 and 1
     1073        else if(levelset[0]==0 && levelset[2]==0){ //front is on point 0 and 1
    10741074                xyz_zero[3*0+0]=xyz_list[2][0];
    10751075                xyz_zero[3*0+1]=xyz_list[2][1];
     
    10811081                xyz_zero[3*1+2]=xyz_list[0][2];
    10821082        }
    1083         else if(levelset[1]==0. && levelset[2]==0.){ //front is on point 0 and 1
     1083        else if(levelset[1]==0 && levelset[2]==0){ //front is on point 0 and 1
    10841084                xyz_zero[3*0+0]=xyz_list[1][0];
    10851085                xyz_zero[3*0+1]=xyz_list[1][1];
     
    14511451
    14521452        /*Go through all the input objects, and find the one corresponding to enum_type, if it exists: */
    1453         if (enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum)
     1453        if (enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyDbarEnum)
    14541454         input=this->material->inputs->GetInput(enum_type);
    14551455        else
     
    15671567                                        break;
    15681568                                case MaterialsRheologyBbarEnum:
    1569                                 case MaterialsRheologyZbarEnum:
     1569                                case MaterialsRheologyDbarEnum:
    15701570                                        /*Material will take care of it*/ break;
    15711571                                default:
     
    17671767void  Tria::InputUpdateFromVector(IssmDouble* vector, int name, int type){
    17681768
     1769
    17691770        /*Check that name is an element input*/
    17701771        if (!IsInput(name)) return;
     
    17801781                }
    17811782                /*update input*/
    1782                 if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheologyZEnum || name==MaterialsRheologyZbarEnum){
     1783                if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheologyDEnum || name==MaterialsRheologyDbarEnum){
    17831784                        material->inputs->AddInput(new TriaInput(name,values,P1Enum));
    17841785                }
     
    17931794                }
    17941795                /*update input*/
    1795                 if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheologyZEnum || name==MaterialsRheologyZbarEnum){
     1796                if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheologyDEnum || name==MaterialsRheologyDbarEnum){
    17961797                        material->inputs->AddInput(new TriaInput(name,values,P1Enum));
    17971798                }
     
    19611962                                name==FrictionCoefficientEnum ||
    19621963                                name==MaterialsRheologyBbarEnum ||
    1963                                 name==MaterialsRheologyZbarEnum ||
     1964                                name==MaterialsRheologyDbarEnum ||
    19641965                                name==GradientEnum ||
    19651966                                name==OldGradientEnum ||
     
    27932794                        *presponse=this->material->GetBbar();
    27942795                        break;
    2795                 case MaterialsRheologyZbarEnum:
    2796                         *presponse=this->material->GetZbar();
    2797                         break;
     2796
    27982797                case VelEnum:{
    27992798
     
    37313730        for(int i=0;i<num_controls;i++){
    37323731
    3733                 if(control_type[i]==MaterialsRheologyBbarEnum || control_type[i]==MaterialsRheologyZbarEnum){
     3732                if(control_type[i]==MaterialsRheologyBbarEnum || control_type[i]==MaterialsRheologyDbarEnum){
    37343733                        input=(Input*)material->inputs->GetInput(control_type[i]); _assert_(input);
    37353734                }
     
    37583757        Input* input=NULL;
    37593758
    3760         if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum){
     3759        if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyDbarEnum){
    37613760                input=(Input*)material->inputs->GetInput(enum_type);
    37623761        }
     
    37763775        Input* input=NULL;
    37773776
    3778         if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum){
     3777        if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyDbarEnum){
    37793778                input=(Input*)material->inputs->GetInput(enum_type);
    37803779        }
     
    37953794        Input* input=NULL;
    37963795
    3797         if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyZbarEnum){
     3796        if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyDbarEnum){
    37983797                input=(Input*)material->inputs->GetInput(enum_type);
    37993798        }
     
    38263825                        GradjBSSA(gradient,control_index);
    38273826                        break;
    3828                 case MaterialsRheologyZbarEnum:
    3829                         GradjZSSA(gradient,control_index);
     3827                case MaterialsRheologyDbarEnum:
     3828                        GradjDSSA(gradient,control_index);
    38303829                        break;
    38313830                case BalancethicknessThickeningRateEnum:
     
    39193918}
    39203919/*}}}*/
    3921 /*FUNCTION Tria::GradjZGradient{{{*/
    3922 void  Tria::GradjZGradient(Vector<IssmDouble>* gradient,int weight_index,int control_index){
    3923 
    3924         int        i;
    3925         int        vertexpidlist[NUMVERTICES];
    3926         IssmDouble Jdet,weight;
    3927         IssmDouble xyz_list[NUMVERTICES][3];
    3928         IssmDouble dbasis[NDOF2][NUMVERTICES];
    3929         IssmDouble dk[NDOF2];
    3930         IssmDouble grade_g[NUMVERTICES]={0.0};
    3931         GaussTria  *gauss=NULL;
    3932 
    3933         /*Retrieve all inputs we will be needing: */
    3934         GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES);
    3935         GradientIndexing(&vertexpidlist[0],control_index);
    3936         Input* rheologyz_input=material->inputs->GetInput(MaterialsRheologyZbarEnum); _assert_(rheologyz_input);
    3937         Input* weights_input=inputs->GetInput(InversionCostFunctionsCoefficientsEnum);                _assert_(weights_input);
    3938 
    3939         /* Start  looping on the number of gaussian points: */
    3940         gauss=new GaussTria(2);
    3941         for(int ig=gauss->begin();ig<gauss->end();ig++){
    3942 
    3943                 gauss->GaussPoint(ig);
    3944 
    3945                 GetJacobianDeterminant(&Jdet, &xyz_list[0][0],gauss);
    3946                 GetNodalFunctionsDerivatives(&dbasis[0][0],&xyz_list[0][0],gauss);
    3947                 weights_input->GetInputValue(&weight,gauss,weight_index);
    3948 
    3949                 /*Build alpha_complement_list: */
    3950                 rheologyz_input->GetInputDerivativeValue(&dk[0],&xyz_list[0][0],gauss);
    3951 
    3952                 /*Build gradje_g_gaussian vector (actually -dJ/ddrag): */
    3953                 for (i=0;i<NUMVERTICES;i++) grade_g[i]+=-weight*Jdet*gauss->weight*(dbasis[0][i]*dk[0]+dbasis[1][i]*dk[1]);
    3954         }
    3955         gradient->SetValues(NUMVERTICES,vertexpidlist,grade_g,ADD_VAL);
    3956 
    3957         /*Clean up and return*/
    3958         delete gauss;
    3959 }
    3960 /*}}}*/
    39613920/*FUNCTION Tria::GradjBSSA{{{*/
    39623921void  Tria::GradjBSSA(Vector<IssmDouble>* gradient,int control_index){
     
    40163975}
    40173976/*}}}*/
    4018 /*FUNCTION Tria::GradjZSSA{{{*/
    4019 void  Tria::GradjZSSA(Vector<IssmDouble>* gradient,int control_index){
     3977/*FUNCTION Tria::GradjDSSA{{{*/
     3978void  Tria::GradjDSSA(Vector<IssmDouble>* gradient,int control_index){
    40203979
    40213980        /*Intermediaries*/
     
    40243983        IssmDouble vx,vy,lambda,mu,thickness,Jdet;
    40253984        IssmDouble viscosity_complement;
    4026         IssmDouble dvx[NDOF2],dvy[NDOF2],dadjx[NDOF2],dadjy[NDOF2],dZ[NDOF2];
     3985        IssmDouble dvx[NDOF2],dvy[NDOF2],dadjx[NDOF2],dadjy[NDOF2];
    40273986        IssmDouble xyz_list[NUMVERTICES][3];
    40283987        IssmDouble basis[3],epsilon[3];
     
    40403999        Input* adjointx_input=inputs->GetInput(AdjointxEnum);                       _assert_(adjointx_input);
    40414000        Input* adjointy_input=inputs->GetInput(AdjointyEnum);                       _assert_(adjointy_input);
    4042         Input* rheologyz_input=material->inputs->GetInput(MaterialsRheologyZbarEnum); _assert_(rheologyz_input);
     4001        Input* rheologyd_input=material->inputs->GetInput(MaterialsRheologyDbarEnum); _assert_(rheologyd_input);
    40434002
    40444003        /* Start  looping on the number of gaussian points: */
     
    40494008
    40504009                thickness_input->GetInputValue(&thickness,gauss);
    4051                 rheologyz_input->GetInputDerivativeValue(&dZ[0],&xyz_list[0][0],gauss);
    40524010                vx_input->GetInputDerivativeValue(&dvx[0],&xyz_list[0][0],gauss);
    40534011                vy_input->GetInputDerivativeValue(&dvy[0],&xyz_list[0][0],gauss);
     
    40564014
    40574015                this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
    4058                 material->GetViscosityZComplement(&viscosity_complement,&epsilon[0]);
     4016                material->GetViscosityDComplement(&viscosity_complement,&epsilon[0]);
    40594017
    40604018                GetJacobianDeterminant(&Jdet, &xyz_list[0][0],gauss);
     
    54215379
    54225380        /*Get input (either in element or material)*/
    5423         if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheologyZbarEnum){
     5381        if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheologyDbarEnum){
    54245382                input=(Input*)material->inputs->GetInput(control_enum); _assert_(input);
    54255383        }
     
    54565414        new_input = new TriaInput(control_enum,values,P1Enum);
    54575415
    5458         if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheologyZbarEnum){
     5416        if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheologyDbarEnum){
    54595417                input=(Input*)material->inputs->GetInput(control_enum); _assert_(input);
    54605418        }
  • issm/trunk-jpl/src/c/classes/Elements/Tria.h

    r16125 r16146  
    139139                void       Gradj(Vector<IssmDouble>* gradient,int control_type,int control_index);
    140140                void       GradjBGradient(Vector<IssmDouble>* gradient,int weight_index,int control_index);
    141                 void       GradjZGradient(Vector<IssmDouble>* gradient,int weight_index,int control_index);
     141                void       GradjDGradient(Vector<IssmDouble>* gradient,int weight_index,int control_index);
    142142                void       GradjBSSA(Vector<IssmDouble>* gradient,int control_index);
    143                 void       GradjZSSA(Vector<IssmDouble>* gradient,int control_index);
     143                void       GradjDSSA(Vector<IssmDouble>* gradient,int control_index);
    144144                void       GradjDragSSA(Vector<IssmDouble>* gradient,int control_index);
    145145                void       GradjDragFS(Vector<IssmDouble>* gradient,int control_index);
  • issm/trunk-jpl/src/c/classes/Materials/Matdamageice.cpp

    r16125 r16146  
    161161}
    162162/*}}}*/
     163/*FUNCTION Matdamageice::GetD {{{*/
     164IssmDouble Matdamageice::GetD(){
     165
     166        /*Output*/
     167        IssmDouble D;
     168
     169        inputs->GetInputAverage(&D,MaterialsRheologyDEnum);
     170        return D;
     171}
     172/*}}}*/
    163173/*FUNCTION Matdamageice::GetZ {{{*/
    164174IssmDouble Matdamageice::GetZ(){
    165175
    166176        /*Output*/
    167         IssmDouble Z;
    168 
    169         inputs->GetInputAverage(&Z,MaterialsRheologyZEnum);
    170         return Z;
     177        IssmDouble D;
     178
     179        inputs->GetInputAverage(&D,MaterialsRheologyDEnum);
     180        return 1/(1-D);
    171181}
    172182/*}}}*/
     
    175185
    176186        /*Output*/
    177         IssmDouble Zbar;
    178 
    179         inputs->GetInputAverage(&Zbar,MaterialsRheologyZbarEnum);
    180         return Zbar;
     187        IssmDouble Dbar;
     188        inputs->GetInputAverage(&Dbar,MaterialsRheologyDbarEnum);
     189        return 1/(1-Dbar);
     190}
     191/*}}}*/
     192/*FUNCTION Matdamageice::GetDbar {{{*/
     193IssmDouble Matdamageice::GetDbar(){
     194
     195        /*Output*/
     196        IssmDouble Dbar;
     197        inputs->GetInputAverage(&Dbar,MaterialsRheologyDbarEnum);
     198        return Dbar;
    181199}
    182200/*}}}*/
     
    216234void  Matdamageice::GetViscosity2d(IssmDouble* pviscosity, IssmDouble* epsilon){
    217235        /*From a string tensor and a material object, return viscosity, using Glen's flow law.
    218                                                                                                    Z * B
     236                                                                                                   (1-D) * B
    219237          viscosity= -------------------------------------------------------------------
    220238                                                  2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
     
    235253        /*Intermediary: */
    236254        IssmDouble A,e;
    237         IssmDouble Btmp,B,n,Z;
     255        IssmDouble Btmp,B,n,D;
    238256
    239257        /*Get B and n*/
    240258        Btmp=GetBbar();
    241         Z=GetZbar();
     259        D=GetDbar();
    242260        n=GetN();
    243         B=Z*Btmp;
     261        B=(1-D)*Btmp;
    244262
    245263        if (n==1){
     
    284302        /*Return viscosity accounting for steady state power law creep [Thomas and SSA, 1982]:
    285303         *
    286          *                                               B
     304         *                                             (1-D)* B
    287305         * viscosity3d= -------------------------------------------------------------------
    288306         *                      2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
     
    303321        /*Intermediaries: */
    304322        IssmDouble A,e;
    305         IssmDouble B,n,Z;
     323        IssmDouble B,n,D;
    306324
    307325        /*Get B, Z and n*/
    308326        n=GetN();
    309         Z=GetZ();
    310         B=Z*GetB();
     327        D=GetD();
     328        B=(1-D)*GetB();
    311329
    312330        if (n==1){
     
    355373        /*Return viscosity accounting for steady state power law creep [Thomas and SSA, 1982]:
    356374         *
    357          *                                          B
     375         *                                             (1-D)* B
    358376         * viscosity3d= -------------------------------------------------------------------
    359377         *                   2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
     
    374392        /*Intermediaries: */
    375393        IssmDouble A,e;
    376         IssmDouble B,n,Z;
     394        IssmDouble B,n,D;
    377395        IssmDouble eps0;
    378396
     
    380398        eps0=pow(10.,(IssmDouble)-27);
    381399        n=GetN();
    382         Z=GetZ();
    383         B=Z*GetB();
     400        D=GetD();
     401        B=(1-D)*GetB();
    384402
    385403        if (n==1){
     
    480498}
    481499/*}}}*/
    482 /*FUNCTION Matdamageice::GetViscosityZComplement {{{*/
    483 void  Matdamageice::GetViscosityZComplement(IssmDouble* pviscosity_complement, IssmDouble* epsilon){
    484         /*Return viscosity derivative for control method d(mu)/dZ:
     500/*FUNCTION Matdamageice::GetViscosityDComplement {{{*/
     501void  Matdamageice::GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* epsilon){
     502        /*Return viscosity derivative for control method d(mu)/dD:
    485503         *
    486504         *                                                                                             B
    487          * dviscosity= -------------------------------------------------------------------
     505         * dviscosity= - -------------------------------------------------------------------
    488506         *                                2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
    489507         *
     
    515533                if(A==0){
    516534                        /*Maximum viscosity_complement for 0 shear areas: */
    517                         viscosity_complement=2.25*pow(10.,17);
     535                        viscosity_complement=- 2.25*pow(10.,17);
    518536                }
    519537                else{
    520538                        e=(n-1)/(2*n);
    521539
    522                         viscosity_complement=B/(2*pow(A,e));
     540                        viscosity_complement=- B/(2*pow(A,e));
    523541                }
    524542        }
    525543        else{
    526                 viscosity_complement=4.5*pow(10.,17);
     544                viscosity_complement=- 4.5*pow(10.,17);
    527545        }
    528546
     
    530548        _assert_(B>0);
    531549        _assert_(n>0);
    532         _assert_(viscosity_complement>0);
     550        _assert_(viscosity_complement<0);
    533551
    534552        /*Return: */
     
    779797                }
    780798
    781                 /*Get Z*/
    782                 if (iomodel->Data(MaterialsRheologyZEnum)) {
    783                         for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheologyZEnum)[iomodel->elements[num_vertices*index+i]-1];
    784                         this->inputs->AddInput(new TriaInput(MaterialsRheologyZbarEnum,nodeinputs,P1Enum));
     799                /*Get D*/
     800                if (iomodel->Data(MaterialsRheologyDEnum)) {
     801                        for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheologyDEnum)[iomodel->elements[num_vertices*index+i]-1];
     802                        this->inputs->AddInput(new TriaInput(MaterialsRheologyDbarEnum,nodeinputs,P1Enum));
    785803                }
    786804
     
    799817                                                }
    800818                                                break;
    801                                         case MaterialsRheologyZbarEnum:
    802                                                 if (iomodel->Data(MaterialsRheologyZEnum)){
    803                                                         _assert_(iomodel->Data(MaterialsRheologyZEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));
    804                                                         for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(MaterialsRheologyZEnum)[iomodel->elements[num_vertices*index+j]-1];
     819                                        case MaterialsRheologyDbarEnum:
     820                                                if (iomodel->Data(MaterialsRheologyDEnum)){
     821                                                        _assert_(iomodel->Data(MaterialsRheologyDEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));
     822                                                        for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(MaterialsRheologyDEnum)[iomodel->elements[num_vertices*index+j]-1];
    805823                                                        for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->Data(InversionMinParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
    806824                                                        for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->Data(InversionMaxParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
    807                                                         this->inputs->AddInput(new ControlInput(MaterialsRheologyZbarEnum,TriaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
     825                                                        this->inputs->AddInput(new ControlInput(MaterialsRheologyDbarEnum,TriaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    808826                                                }
    809827                                                break;
     
    837855                }
    838856
    839                 /*Get Z*/
    840                 if (iomodel->Data(MaterialsRheologyZEnum)) {
    841                         for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheologyZEnum)[iomodel->elements[num_vertices*index+i]-1];
    842                         this->inputs->AddInput(new PentaInput(MaterialsRheologyZEnum,nodeinputs,P1Enum));
     857                /*Get D*/
     858                if (iomodel->Data(MaterialsRheologyDEnum)) {
     859                        for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheologyDEnum)[iomodel->elements[num_vertices*index+i]-1];
     860                        this->inputs->AddInput(new PentaInput(MaterialsRheologyDEnum,nodeinputs,P1Enum));
    843861                }
    844862
     
    857875                                                }
    858876                                                break;
    859                                         case MaterialsRheologyZbarEnum:
    860                                                 if (iomodel->Data(MaterialsRheologyZEnum)){
    861                                                         _assert_(iomodel->Data(MaterialsRheologyZEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));
    862                                                         for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(MaterialsRheologyZEnum)[iomodel->elements[num_vertices*index+j]-1];
     877                                        case MaterialsRheologyDbarEnum:
     878                                                if (iomodel->Data(MaterialsRheologyDEnum)){
     879                                                        _assert_(iomodel->Data(MaterialsRheologyDEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));
     880                                                        for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(MaterialsRheologyDEnum)[iomodel->elements[num_vertices*index+j]-1];
    863881                                                        for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->Data(InversionMinParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
    864882                                                        for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->Data(InversionMaxParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
    865                                                         this->inputs->AddInput(new ControlInput(MaterialsRheologyZEnum,PentaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
     883                                                        this->inputs->AddInput(new ControlInput(MaterialsRheologyDEnum,PentaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
    866884                                                }
    867885                                                break;
     
    885903                                name==MaterialsRheologyBbarEnum ||
    886904                                name==MaterialsRheologyNEnum ||
    887                                 name==MaterialsRheologyZEnum ||
    888                                 name==MaterialsRheologyZbarEnum
     905                                name==MaterialsRheologyDEnum ||
     906                                name==MaterialsRheologyDbarEnum
    889907                ){
    890908                return true;
  • issm/trunk-jpl/src/c/classes/Materials/Matdamageice.h

    r16125 r16146  
    5050                void   GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon);
    5151                void   GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
    52                 void   GetViscosityZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
     52                void   GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
    5353                void   GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
    5454                void   GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon);
     
    5959                IssmDouble GetN();
    6060                IssmDouble GetZ();
     61                IssmDouble GetD();
    6162                IssmDouble GetZbar();
     63                IssmDouble GetDbar();
    6264                bool   IsInput(int name);
    6365                /*}}}*/
  • issm/trunk-jpl/src/c/classes/Materials/Material.h

    r15654 r16146  
    3030                virtual void       GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon)=0;
    3131                virtual void       GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon)=0;
    32                 virtual void       GetViscosityZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon)=0;
     32                virtual void       GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon)=0;
    3333                virtual void       GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon)=0;
    3434                virtual void       GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon)=0;
     
    3939                virtual IssmDouble GetN()=0;
    4040                virtual IssmDouble GetZ()=0;
     41                virtual IssmDouble GetD()=0;
    4142                virtual IssmDouble GetZbar()=0;
     43                virtual IssmDouble GetDbar()=0;
    4244
    4345};
  • issm/trunk-jpl/src/c/classes/Materials/Matice.h

    r16125 r16146  
    5757                void       GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon);
    5858                void       GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
    59                 void GetViscosityZComplement(IssmDouble*, IssmDouble*){_error_("not supported");};
     59                void       GetViscosityDComplement(IssmDouble*, IssmDouble*){_error_("not supported");};
    6060                void       GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
    6161                void       GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon);
     
    6565                IssmDouble GetBbar();
    6666                IssmDouble GetZ(){_error_("not supported");};
     67                IssmDouble GetD(){_error_("not supported");};
    6768                IssmDouble GetZbar(){_error_("not supported");};
     69                IssmDouble GetDbar(){_error_("not supported");};
    6870                IssmDouble GetN();
    6971                bool       IsInput(int name);
  • issm/trunk-jpl/src/c/classes/Materials/Matpar.h

    r16125 r16146  
    8787                void       GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon){_error_("not supported");};
    8888                void       GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon){_error_("not supported");};
    89                 void       GetViscosityZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon){_error_("not supported");};
     89                void       GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon){_error_("not supported");};
    9090                void       GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon){_error_("not supported");};
    9191                void       GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon){_error_("not supported");};
     
    9595                IssmDouble GetBbar(){_error_("not supported");};
    9696                IssmDouble GetN(){_error_("not supported");};
     97                IssmDouble GetD(){_error_("not supported");};
    9798                IssmDouble GetZ(){_error_("not supported");};
    9899                IssmDouble GetZbar(){_error_("not supported");};
     100                IssmDouble GetDbar(){_error_("not supported");};
    99101                /*}}}*/
    100102                /*Numerics: {{{*/
  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Control/UpdateElementsAndMaterialsControl.cpp

    r16050 r16146  
    4040                        case FrictionCoefficientEnum:   iomodel->FetchData(1,FrictionCoefficientEnum); break;
    4141                        case MaterialsRheologyBbarEnum: iomodel->FetchData(1,MaterialsRheologyBEnum); break;
    42                         case MaterialsRheologyZbarEnum: iomodel->FetchData(1,MaterialsRheologyZEnum); break;
     42                        case MaterialsRheologyDbarEnum: iomodel->FetchData(1,MaterialsRheologyDEnum); break;
    4343                        default: _error_("Control " << EnumToStringx(reCast<int,IssmDouble>(iomodel->Data(InversionControlParametersEnum)[i])) << " not implemented yet");
    4444                }
     
    5959
    6060        /*Free data: */
    61         iomodel->DeleteData(4+7,InversionControlParametersEnum,InversionCostFunctionsCoefficientsEnum,InversionMinParametersEnum,InversionMaxParametersEnum,BalancethicknessThickeningRateEnum,VxEnum,VyEnum,ThicknessEnum,FrictionCoefficientEnum,MaterialsRheologyBEnum,MaterialsRheologyZEnum);
     61        iomodel->DeleteData(4+7,InversionControlParametersEnum,InversionCostFunctionsCoefficientsEnum,InversionMinParametersEnum,InversionMaxParametersEnum,BalancethicknessThickeningRateEnum,VxEnum,VyEnum,ThicknessEnum,FrictionCoefficientEnum,MaterialsRheologyBEnum,MaterialsRheologyDEnum);
    6262}
  • issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateElementsVerticesAndMaterials.cpp

    r15941 r16146  
    6464                        break;
    6565                case MatdamageiceEnum:
    66                         iomodel->FetchData(3,MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheologyZEnum);
     66                        iomodel->FetchData(3,MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheologyDEnum);
    6767                        for (i=0;i<iomodel->numberofelements;i++) if(iomodel->my_elements[i]) materials->AddObject(new Matdamageice(i+1,i,iomodel));
    6868                        break;
     
    7373        /*Free data: */
    7474        iomodel->DeleteData(8,MeshUpperelementsEnum,MeshLowerelementsEnum,
    75                                 MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheologyZEnum,InversionControlParametersEnum,InversionMinParametersEnum,
     75                                MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheologyDEnum,InversionControlParametersEnum,InversionMinParametersEnum,
    7676                                InversionMaxParametersEnum);
    7777
  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Stressbalance/UpdateElementsStressbalance.cpp

    r16144 r16146  
    102102        iomodel->FetchDataToInput(elements,LoadingforceYEnum);
    103103        if(materials_type==MatdamageiceEnum){
    104                 iomodel->FetchDataToInput(elements,MaterialsRheologyZEnum);
     104                iomodel->FetchDataToInput(elements,MaterialsRheologyDEnum);
    105105        }
    106106        if(iomodel->dim==3){
  • issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

    r16142 r16146  
    154154        MaterialsRheologyZEnum,
    155155        MaterialsRheologyZbarEnum,
     156        MaterialsRheologyDEnum,
     157        MaterialsRheologyDbarEnum,
    156158        MaterialsRhoIceEnum,
    157159        MaterialsRhoWaterEnum,
  • issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

    r16142 r16146  
    162162                case MaterialsRheologyZEnum : return "MaterialsRheologyZ";
    163163                case MaterialsRheologyZbarEnum : return "MaterialsRheologyZbar";
     164                case MaterialsRheologyDEnum : return "MaterialsRheologyD";
     165                case MaterialsRheologyDbarEnum : return "MaterialsRheologyDbar";
    164166                case MaterialsRhoIceEnum : return "MaterialsRhoIce";
    165167                case MaterialsRhoWaterEnum : return "MaterialsRhoWater";
  • issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

    r16142 r16146  
    165165              else if (strcmp(name,"MaterialsRheologyZ")==0) return MaterialsRheologyZEnum;
    166166              else if (strcmp(name,"MaterialsRheologyZbar")==0) return MaterialsRheologyZbarEnum;
     167              else if (strcmp(name,"MaterialsRheologyD")==0) return MaterialsRheologyDEnum;
     168              else if (strcmp(name,"MaterialsRheologyDbar")==0) return MaterialsRheologyDbarEnum;
    167169              else if (strcmp(name,"MaterialsRhoIce")==0) return MaterialsRhoIceEnum;
    168170              else if (strcmp(name,"MaterialsRhoWater")==0) return MaterialsRhoWaterEnum;
     
    258260              else if (strcmp(name,"GiaLithosphereThickness")==0) return GiaLithosphereThicknessEnum;
    259261              else if (strcmp(name,"Thickness")==0) return ThicknessEnum;
    260               else if (strcmp(name,"TimesteppingStartTime")==0) return TimesteppingStartTimeEnum;
    261               else if (strcmp(name,"TimesteppingFinalTime")==0) return TimesteppingFinalTimeEnum;
    262262         else stage=3;
    263263   }
    264264   if(stage==3){
    265               if (strcmp(name,"TimesteppingCflCoefficient")==0) return TimesteppingCflCoefficientEnum;
     265              if (strcmp(name,"TimesteppingStartTime")==0) return TimesteppingStartTimeEnum;
     266              else if (strcmp(name,"TimesteppingFinalTime")==0) return TimesteppingFinalTimeEnum;
     267              else if (strcmp(name,"TimesteppingCflCoefficient")==0) return TimesteppingCflCoefficientEnum;
    266268              else if (strcmp(name,"TimesteppingTimeAdapt")==0) return TimesteppingTimeAdaptEnum;
    267269              else if (strcmp(name,"TimesteppingTimeStep")==0) return TimesteppingTimeStepEnum;
     
    381383              else if (strcmp(name,"Penpair")==0) return PenpairEnum;
    382384              else if (strcmp(name,"Penta")==0) return PentaEnum;
    383               else if (strcmp(name,"PentaInput")==0) return PentaInputEnum;
    384               else if (strcmp(name,"Profiler")==0) return ProfilerEnum;
    385385         else stage=4;
    386386   }
    387387   if(stage==4){
    388               if (strcmp(name,"MatrixParam")==0) return MatrixParamEnum;
     388              if (strcmp(name,"PentaInput")==0) return PentaInputEnum;
     389              else if (strcmp(name,"Profiler")==0) return ProfilerEnum;
     390              else if (strcmp(name,"MatrixParam")==0) return MatrixParamEnum;
    389391              else if (strcmp(name,"NodeSId")==0) return NodeSIdEnum;
    390392              else if (strcmp(name,"VectorParam")==0) return VectorParamEnum;
     
    504506              else if (strcmp(name,"P0")==0) return P0Enum;
    505507              else if (strcmp(name,"P1")==0) return P1Enum;
    506               else if (strcmp(name,"P1DG")==0) return P1DGEnum;
    507               else if (strcmp(name,"P1bubble")==0) return P1bubbleEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"P1bubblecondensed")==0) return P1bubblecondensedEnum;
     511              if (strcmp(name,"P1DG")==0) return P1DGEnum;
     512              else if (strcmp(name,"P1bubble")==0) return P1bubbleEnum;
     513              else if (strcmp(name,"P1bubblecondensed")==0) return P1bubblecondensedEnum;
    512514              else if (strcmp(name,"P2")==0) return P2Enum;
    513515              else if (strcmp(name,"P2xP1")==0) return P2xP1Enum;
  • issm/trunk-jpl/src/m/classes/inversion.m

    r15860 r16146  
    116116                        md = checkfield(md,'inversion.incomplete_adjoint','values',[0 1]);
    117117                        md = checkfield(md,'inversion.control_parameters','cell',1,'values',...
    118                                 {'BalancethicknessThickeningRate' 'FrictionCoefficient' 'MaterialsRheologyBbar' 'MaterialsRheologyZbar' 'Vx' 'Vy' 'Thickness'});
     118                                {'BalancethicknessThickeningRate' 'FrictionCoefficient' 'MaterialsRheologyBbar' 'MaterialsRheologyDbar' 'Vx' 'Vy' 'Thickness'});
    119119                        md = checkfield(md,'inversion.nsteps','numel',1,'>=',1);
    120120                        md = checkfield(md,'inversion.maxiter_per_step','size',[md.inversion.nsteps 1],'>=',0);
  • issm/trunk-jpl/src/m/classes/inversion.py

    r15860 r16146  
    124124                md = checkfield(md,'inversion.tao','values',[0,1])
    125125                md = checkfield(md,'inversion.incomplete_adjoint','values',[0,1])
    126                 md = checkfield(md,'inversion.control_parameters','cell',1,'values',['BalancethicknessThickeningRate','FrictionCoefficient','MaterialsRheologyBbar','MaterialsRheologyZbar','Vx','Vy'])
     126                md = checkfield(md,'inversion.control_parameters','cell',1,'values',['BalancethicknessThickeningRate','FrictionCoefficient','MaterialsRheologyBbar','MaterialsRheologyDbar','Vx','Vy'])
    127127                md = checkfield(md,'inversion.nsteps','numel',[1],'>=',1)
    128128                md = checkfield(md,'inversion.maxiter_per_step','size',[md.inversion.nsteps],'>=',0)
  • issm/trunk-jpl/src/m/classes/matdamageice.m

    r15825 r16146  
    1919                rheology_B   = NaN;
    2020                rheology_n   = NaN;
    21                 rheology_Z   = NaN;
     21                rheology_D   = NaN;
    2222                rheology_law = '';
    2323
     
    101101                        md = checkfield(md,'materials.rheology_B','>',0,'size',[md.mesh.numberofvertices 1]);
    102102                        md = checkfield(md,'materials.rheology_n','>',0,'size',[md.mesh.numberofelements 1]);
    103                         md = checkfield(md,'materials.rheology_Z','>',0,'size',[md.mesh.numberofvertices 1]);
     103                        md = checkfield(md,'materials.rheology_D','>',0,'size',[md.mesh.numberofvertices 1]);
    104104                        md = checkfield(md,'materials.rheology_law','values',{'None' 'Paterson' 'Arrhenius' 'LliboutryDuval'});
    105105
     
    128128                        fielddisplay(obj,'rheology_B','flow law parameter [Pa/s^(1/n)]');
    129129                        fielddisplay(obj,'rheology_n','Glen''s flow law exponent');
    130                         fielddisplay(obj,'rheology_Z','rheology multiplier');
     130                        fielddisplay(obj,'rheology_D','damage tensor (scalar)');
    131131                        fielddisplay(obj,'rheology_law',['law for the temperature dependance of the rheology: ''None'', ''Paterson'', ''Arrhenius'' or ''LliboutryDuval''']);
    132132                        fielddisplay(obj,'lithosphere_shear_modulus','Lithosphere shear modulus [Pa]');
     
    151151                        WriteData(fid,'object',obj,'class','materials','fieldname','rheology_B','format','DoubleMat','mattype',1);
    152152                        WriteData(fid,'object',obj,'class','materials','fieldname','rheology_n','format','DoubleMat','mattype',2);
    153                         WriteData(fid,'object',obj,'class','materials','fieldname','rheology_Z','format','DoubleMat','mattype',1);
     153                        WriteData(fid,'object',obj,'class','materials','fieldname','rheology_D','format','DoubleMat','mattype',1);
    154154                        WriteData(fid,'data',StringToEnum(obj.rheology_law),'enum',MaterialsRheologyLawEnum(),'format','Integer');
    155155
  • issm/trunk-jpl/src/m/classes/matdamageice.py

    r15852 r16146  
    2727                self.rheology_B   = float('NaN')
    2828                self.rheology_n   = float('NaN')
    29                 self.rheology_Z   = float('NaN')
     29                self.rheology_D   = float('NaN')
    3030                self.rheology_law = ''
    3131
     
    5858                s+="%s\n" % fielddisplay(self,"rheology_B","flow law parameter [Pa/s^(1/n)]")
    5959                s+="%s\n" % fielddisplay(self,"rheology_n","Glen's flow law exponent")
    60                 s+="%s\n" % fielddisplay(self,"rheology_Z","rheology multiplier")
     60                s+="%s\n" % fielddisplay(self,"rheology_D","damage tensor (scalar for now)")
    6161                s+="%s\n" % fielddisplay(self,"rheology_law","law for the temperature dependance of the rheology: 'None', 'Paterson', 'Arrhenius' or 'LliboutryDuval'")
    6262                s+="%s\n" % fielddisplay(self,"lithosphere_shear_modulus","Lithosphere shear modulus [Pa]")
     
    119119                md = checkfield(md,'materials.rheology_B','>',0,'size',[md.mesh.numberofvertices])
    120120                md = checkfield(md,'materials.rheology_n','>',0,'size',[md.mesh.numberofelements])
    121                 md = checkfield(md,'materials.rheology_Z','>',0,'size',[md.mesh.numberofvertices])
     121                md = checkfield(md,'materials.rheology_D','>',0,'size',[md.mesh.numberofvertices])
    122122                md = checkfield(md,'materials.rheology_law','values',['None','Paterson','Arrhenius','LliboutryDuval'])
    123123                md = checkfield(md,'materials.lithosphere_shear_modulus','>',0,'numel',[1]);
     
    143143                WriteData(fid,'object',self,'class','materials','fieldname','rheology_B','format','DoubleMat','mattype',1)
    144144                WriteData(fid,'object',self,'class','materials','fieldname','rheology_n','format','DoubleMat','mattype',2)
    145                 WriteData(fid,'object',self,'class','materials','fieldname','rheology_Z','format','DoubleMat','mattype',1)
     145                WriteData(fid,'object',self,'class','materials','fieldname','rheology_D','format','DoubleMat','mattype',1)
    146146                WriteData(fid,'data',StringToEnum(self.rheology_law)[0],'enum',MaterialsRheologyLawEnum(),'format','Integer')
    147147
  • issm/trunk-jpl/src/m/enum/EnumDefinitions.py

    r16142 r16146  
    154154def MaterialsRheologyZEnum(): return StringToEnum("MaterialsRheologyZ")[0]
    155155def MaterialsRheologyZbarEnum(): return StringToEnum("MaterialsRheologyZbar")[0]
     156def MaterialsRheologyDEnum(): return StringToEnum("MaterialsRheologyD")[0]
     157def MaterialsRheologyDbarEnum(): return StringToEnum("MaterialsRheologyDbar")[0]
    156158def MaterialsRhoIceEnum(): return StringToEnum("MaterialsRhoIce")[0]
    157159def MaterialsRhoWaterEnum(): return StringToEnum("MaterialsRhoWater")[0]
  • issm/trunk-jpl/test/NightlyRun/test270.m

    r15771 r16146  
    55md.materials.rheology_B=paterson(md.initialization.temperature);
    66md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    7 md.materials.rheology_Z=0.5*ones(md.mesh.numberofvertices,1);
     7md.materials.rheology_D=0.5*ones(md.mesh.numberofvertices,1);
    88md=setflowequation(md,'SSA','all');
    99md.cluster=generic('name',oshostname(),'np',3);
  • issm/trunk-jpl/test/NightlyRun/test270.py

    r15771 r16146  
    1717md.materials.rheology_B=paterson(md.initialization.temperature)
    1818md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    19 md.materials.rheology_Z=0.5*numpy.ones((md.mesh.numberofvertices,1))
     19md.materials.rheology_D=0.5*numpy.ones((md.mesh.numberofvertices,1))
    2020md=setflowequation(md,'SSA','all')
    2121md.cluster=generic('name',oshostname(),'np',3)
  • issm/trunk-jpl/test/NightlyRun/test272.m

    r15771 r16146  
    55md.materials.rheology_B=paterson(md.initialization.temperature);
    66md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    7 md.materials.rheology_Z=0.5*ones(md.mesh.numberofvertices,1);
     7md.materials.rheology_D=0.5*ones(md.mesh.numberofvertices,1);
    88md=setflowequation(md,'SSA','all');
     9md.verbose=verbose('control',true);
    910
    1011%control parameters
    1112md.inversion.iscontrol=1;
    12 md.inversion.control_parameters={'MaterialsRheologyZbar'};
    13 md.inversion.min_parameters=10^-13*ones(md.mesh.numberofvertices,1);
    14 md.inversion.max_parameters=ones(md.mesh.numberofvertices,1);
     13md.inversion.control_parameters={'MaterialsRheologyDbar'};
     14md.inversion.min_parameters=zeros(md.mesh.numberofvertices,1);
     15md.inversion.max_parameters=(1-10^-13)*ones(md.mesh.numberofvertices,1);
    1516md.inversion.nsteps=2;
    1617md.inversion.cost_functions=101.*ones(md.inversion.nsteps,1);
     
    2627
    2728%Fields and tolerances to track changes
    28 field_names     ={'Gradient','Misfits','MaterialsRheologyZbar','Pressure','Vel','Vx','Vy'};
     29field_names     ={'Gradient','Misfits','MaterialsRheologyDbar','Pressure','Vel','Vx','Vy'};
    2930field_tolerances={1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12};
    3031field_values={...
    3132   (md.results.StressbalanceSolution.Gradient1),...
    3233   (md.results.StressbalanceSolution.J),...
    33    (md.results.StressbalanceSolution.MaterialsRheologyZbar),...
     34   (md.results.StressbalanceSolution.MaterialsRheologyDbar),...
    3435   (md.results.StressbalanceSolution.Pressure),...
    3536   (md.results.StressbalanceSolution.Vel),...
  • issm/trunk-jpl/test/NightlyRun/test272.py

    r15771 r16146  
    1717md.materials.rheology_B=paterson(md.initialization.temperature)
    1818md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    19 md.materials.rheology_Z=0.5*numpy.ones((md.mesh.numberofvertices,1))
     19md.materials.rheology_D=0.5*numpy.ones((md.mesh.numberofvertices,1))
    2020md=setflowequation(md,'SSA','all')
    2121
    2222#control parameters
    2323md.inversion.iscontrol=1
    24 md.inversion.control_parameters=['MaterialsRheologyZbar']
     24md.inversion.control_parameters=['MaterialsRheologyDbar']
    2525md.inversion.min_parameters=10**-13*numpy.ones((md.mesh.numberofvertices,1))
    2626md.inversion.max_parameters=numpy.ones((md.mesh.numberofvertices,1))
     
    3838
    3939#Fields and tolerances to track changes
    40 field_names     =['Gradient','Misfits','MaterialsRheologyZbar','Pressure','Vel','Vx','Vy']
     40field_names     =['Gradient','Misfits','MaterialsRheologyDbar','Pressure','Vel','Vx','Vy']
    4141field_tolerances=[1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12]
    4242field_values=[\
    4343   md.results.StressbalanceSolution.Gradient1,\
    4444   md.results.StressbalanceSolution.J,\
    45    md.results.StressbalanceSolution.MaterialsRheologyZbar,\
     45   md.results.StressbalanceSolution.MaterialsRheologyDbar,\
    4646   md.results.StressbalanceSolution.Pressure,\
    4747   md.results.StressbalanceSolution.Vel,\
  • issm/trunk-jpl/test/NightlyRun/test274.m

    r16075 r16146  
    66md.materials.rheology_B=paterson(md.initialization.temperature);
    77md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    8 md.materials.rheology_Z=0.5*ones(md.mesh.numberofvertices,1);
     8md.materials.rheology_D=0.5*ones(md.mesh.numberofvertices,1);
    99md=setflowequation(md,'SSA','all');
    1010
  • issm/trunk-jpl/test/NightlyRun/test274.py

    r16075 r16146  
    1919md.materials.rheology_B=paterson(md.initialization.temperature)
    2020md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    21 md.materials.rheology_Z=0.5*numpy.ones((md.mesh.numberofvertices,1))
     21md.materials.rheology_D=0.5*numpy.ones((md.mesh.numberofvertices,1))
    2222md=setflowequation(md,'SSA','all')
    2323
Note: See TracChangeset for help on using the changeset viewer.