Changeset 16167


Ignore:
Timestamp:
09/18/13 12:09:59 (12 years ago)
Author:
Eric.Larour
Message:

CHG: switched off Matdamageice. Damange D and rheology B are now equally treated in the
material Matice. Adjusted all control methods and parameter files accordingly.

Location:
issm/trunk-jpl
Files:
7 deleted
62 edited

Legend:

Unmodified
Added
Removed

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

    r16126 r16167  
    100100                                        ./classes/Materials/Matice.h\
    101101                                        ./classes/Materials/Matice.cpp\
    102                                         ./classes/Materials/Matdamageice.h\
    103                                         ./classes/Materials/Matdamageice.cpp\
    104102                                        ./classes/Materials/Matpar.h\
    105103                                        ./classes/Materials/Matpar.cpp\

  • issm/trunk-jpl/src/c/classes/Elements/Penta.cpp

    r16166 r16167  
    21422142
    21432143        /*Go through all the input objects, and find the one corresponding to enum_type, if it exists: */
    2144         if      (enum_type==MaterialsRheologyBbarEnum)
    2145          input=this->material->inputs->GetInput(MaterialsRheologyBEnum);
    2146         else if (enum_type==MaterialsRheologyZbarEnum)
    2147          input=this->material->inputs->GetInput(MaterialsRheologyZEnum);
    2148         else
    2149          input=this->inputs->GetInput(enum_type);
     2144        if (enum_type==MaterialsRheologyBbarEnum) input=this->material->inputs->GetInput(MaterialsRheologyBEnum);
     2145        else if (enum_type==DamageDbarEnum) input=this->material->inputs->GetInput(DamageDEnum);
     2146        else input=this->inputs->GetInput(enum_type);
    21502147        //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
    21512148        //try and output whatever we can instead of just failing.
    … …  
    22582255                                        }
    22592256                                        break;
    2260                                 case MaterialsRheologyBbarEnum:
    2261                                 case MaterialsRheologyZbarEnum:
    2262                                         /*Material will take care of it*/ break;
     2257                                /*Material will take care of it*/
     2258                                case MaterialsRheologyBbarEnum: break;
     2259                                case DamageDbarEnum:break;
    22632260                                default:
    22642261                                        _error_("Control " << EnumToStringx(reCast<int,IssmDouble>(iomodel->Data(InversionControlParametersEnum)[i])) << " not implemented yet");
    … …  
    26672664                        }
    26682665                        /*update input*/
    2669                         if (name==MaterialsRheologyZEnum || name==MaterialsRheologyZbarEnum){
    2670                                 material->inputs->AddInput(new PentaInput(name,values,P1Enum));
    2671                         }
    2672                         else{
    2673                                 this->inputs->AddInput(new PentaInput(name,values,P1Enum));
    2674                         }
     2666                        this->inputs->AddInput(new PentaInput(name,values,P1Enum));
    26752667                        return;
    26762668
    … …  
    26802672                        }
    26812673                        /*update input*/
    2682                         if (name==MaterialsRheologyZEnum || name==MaterialsRheologyZbarEnum){
    2683                                 material->inputs->AddInput(new PentaInput(name,values,P1Enum));
    2684                         }
    2685                         else{
    2686                                 this->inputs->AddInput(new PentaInput(name,values,P1Enum));
    2687                         }
     2674                        this->inputs->AddInput(new PentaInput(name,values,P1Enum));
    26882675                        return;
    26892676
    … …  
    39803967                        *presponse=this->material->GetBbar();
    39813968                        break;
    3982                 case MaterialsRheologyZbarEnum:
    3983                         *presponse=this->material->GetZbar();
     3969                case DamageDbarEnum:
     3970                        *presponse=this->material->GetDbar();
    39843971                        break;
    39853972                case VelEnum:
    … …  
    54335420                input=(Input*)material->inputs->GetInput(MaterialsRheologyBEnum);
    54345421        }
    5435         else if(enum_type==MaterialsRheologyZbarEnum){
     5422        else if(enum_type==DamageDbarEnum){
    54365423                if(!IsOnBed()) return;
    5437                 input=(Input*)material->inputs->GetInput(MaterialsRheologyZEnum);
    5438         }
    5439 
     5424                input=(Input*)material->inputs->GetInput(DamageDEnum);
     5425        }
    54405426        else{
    54415427                input=inputs->GetInput(enum_type);
    … …  
    54565442                input=(Input*)material->inputs->GetInput(MaterialsRheologyBEnum);
    54575443        }
    5458         else if(enum_type==MaterialsRheologyZbarEnum){
    5459                 input=(Input*)material->inputs->GetInput(MaterialsRheologyZEnum);
     5444        else if(enum_type==DamageDbarEnum){
     5445                input=(Input*)material->inputs->GetInput(DamageDEnum);
    54605446        }
    54615447        else{
    … …  
    54785464                input=(Input*)material->inputs->GetInput(MaterialsRheologyBEnum);
    54795465        }
    5480         else if(enum_type==MaterialsRheologyZbarEnum){
    5481                 input=(Input*)material->inputs->GetInput(MaterialsRheologyZEnum);
     5466        else if(enum_type==DamageDbarEnum){
     5467                input=(Input*)material->inputs->GetInput(DamageDEnum);
    54825468        }
    54835469        else{
    … …  
    55255511                case MaticeEnum:
    55265512                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    5527                         break;
    5528                 case MatdamageiceEnum:
    5529                         this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    5530                         this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     5513                        this->InputDepthAverageAtBase(DamageDEnum,DamageDbarEnum,MaterialsEnum);
    55315514                        break;
    55325515                default:
    … …  
    55415524        /*Delete averaged fields*/
    55425525        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    5543         this->material->inputs->DeleteInput(MaterialsRheologyZbarEnum);
     5526        this->material->inputs->DeleteInput(DamageDbarEnum);
    55445527
    55455528        /*clean up and return*/
    … …  
    61376120        /*If on water, skip grad (=0): */
    61386121        if(NoIceInElement())return;
    6139 
     6122                                       
    61406123        /*First deal with ∂/∂alpha(KU-F)*/
    61416124        switch(control_type){
    … …  
    64046387        /*Depth Average B*/
    64056388        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     6389        this->InputDepthAverageAtBase(DamageDEnum,DamageDbarEnum,MaterialsEnum);
    64066390
    64076391        /*Collapse element to the base*/
    … …  
    64126396        /*delete Average B*/
    64136397        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     6398        this->material->inputs->DeleteInput(DamageDbarEnum);
    64146399
    64156400} /*}}}*/
    … …  
    64206405        if (!IsOnBed()) return;
    64216406
    6422         /*Depth Average B*/
     6407        /*Depth Average B and D*/
    64236408        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     6409        this->InputDepthAverageAtBase(DamageDEnum,DamageDbarEnum,MaterialsEnum);
    64246410
    64256411        /*Collapse element to the base*/
    … …  
    64306416        /*delete Average B*/
    64316417        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     6418        this->material->inputs->DeleteInput(DamageDbarEnum);
    64326419} /*}}}*/
    64336420/*FUNCTION Penta::GradjBbarFS {{{*/
    … …  
    64376424        if (!IsOnBed()) return;
    64386425
    6439         /*Depth Average B*/
     6426        /*Depth Average B and D*/
    64406427        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
     6428        this->InputDepthAverageAtBase(DamageDEnum,DamageDbarEnum,MaterialsEnum);
    64416429
    64426430        /*Collapse element to the base*/
    … …  
    64476435        /*delete Average B*/
    64486436        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
     6437        this->material->inputs->DeleteInput(DamageDbarEnum);
    64496438} /*}}}*/
    64506439/*FUNCTION Penta::InputControlUpdate{{{*/
    … …  
    64666455                        input=(Input*)material->inputs->GetInput(MaterialsRheologyBEnum); _assert_(input);
    64676456                }
    6468                 else if(control_type[i]==MaterialsRheologyZbarEnum){
     6457                else if(control_type[i]==DamageDbarEnum){
    64696458                        if (!IsOnBed()) goto cleanup_and_return;
    6470                         input=(Input*)material->inputs->GetInput(MaterialsRheologyZEnum); _assert_(input);
     6459                        input=(Input*)material->inputs->GetInput(DamageDEnum); _assert_(input);
    64716460                }
    64726461                else{
    … …  
    64836472                        this->InputExtrude(MaterialsRheologyBEnum,MaterialsEnum);
    64846473                }
    6485                 else if(control_type[i]==MaterialsRheologyZbarEnum){
    6486                         this->InputExtrude(MaterialsRheologyZEnum,MaterialsEnum);
     6474                else if(control_type[i]==DamageDbarEnum){
     6475                        this->InputExtrude(DamageDEnum,MaterialsEnum);
    64876476                }
    64886477        }
    … …  
    76957684                case MaticeEnum:
    76967685                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    7697                         break;
    7698                 case MatdamageiceEnum:
    7699                         this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    7700                         this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     7686                        this->InputDepthAverageAtBase(DamageDEnum,DamageDbarEnum,MaterialsEnum);
    77017687                        break;
    77027688                default:
    … …  
    77117697        /*Delete averaged fields*/
    77127698        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    7713         this->material->inputs->DeleteInput(MaterialsRheologyZbarEnum);
     7699        this->material->inputs->DeleteInput(DamageDbarEnum);
    77147700
    77157701        /*clean up and return*/
    … …  
    97379723                case MaticeEnum:
    97389724                        this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    9739                         break;
    9740                 case MatdamageiceEnum:
    9741                         this->InputDepthAverageAtBase(MaterialsRheologyBEnum,MaterialsRheologyBbarEnum,MaterialsEnum);
    9742                         this->InputDepthAverageAtBase(MaterialsRheologyZEnum,MaterialsRheologyZbarEnum,MaterialsEnum);
     9725                        this->InputDepthAverageAtBase(DamageDEnum,DamageDbarEnum,MaterialsEnum);
    97439726                        break;
    97449727                default:
    … …  
    97539736        /*Delete averaged inputs*/
    97549737        this->material->inputs->DeleteInput(MaterialsRheologyBbarEnum);
    9755         this->material->inputs->DeleteInput(MaterialsRheologyZbarEnum);
     9738        this->material->inputs->DeleteInput(DamageDbarEnum);
    97569739
    97579740        /*clean up and return*/

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

    r16164 r16167  
    15611561                                        }
    15621562                                        break;
    1563                                 case MaterialsRheologyBbarEnum:
    1564                                 case DamageDbarEnum:
    1565                                         /*Material will take care of it*/ break;
     1563                                /*Material will take care of it*/
     1564                                case MaterialsRheologyBbarEnum:break;
     1565                                case DamageDbarEnum: break;
    15661566                                default:
    15671567                                        _error_("Control " << EnumToStringx(reCast<int,IssmDouble>(iomodel->Data(InversionControlParametersEnum)[i])) << " not implemented yet");

  • issm/trunk-jpl/src/c/classes/Inputs/ControlInput.cpp

    r15737 r16167  
    6161                /*FUNCTION ControlInput::Echo {{{*/
    6262void ControlInput::Echo(void){
     63        printf("I'm here\n");
    6364        this->DeepEcho();
    6465}

  • issm/trunk-jpl/src/c/classes/Materials/Matice.cpp

    r16125 r16167  
    175175}
    176176/*}}}*/
     177/*FUNCTION Matice::GetD {{{*/
     178IssmDouble Matice::GetD(){
     179
     180        /*Output*/
     181        IssmDouble D;
     182
     183        inputs->GetInputAverage(&D,DamageDEnum);
     184        return D;
     185}
     186/*}}}*/
     187/*FUNCTION Matice::GetZ {{{*/
     188IssmDouble Matice::GetZ(){
     189
     190        /*Output*/
     191        IssmDouble D;
     192
     193        inputs->GetInputAverage(&D,DamageDEnum);
     194        return 1/(1-D);
     195}
     196/*}}}*/
     197/*FUNCTION Matice::GetZbar {{{*/
     198IssmDouble Matice::GetZbar(){
     199
     200        /*Output*/
     201        IssmDouble Dbar;
     202        inputs->GetInputAverage(&Dbar,DamageDbarEnum);
     203        return 1/(1-Dbar);
     204}
     205/*}}}*/
     206/*FUNCTION Matice::GetDbar {{{*/
     207IssmDouble Matice::GetDbar(){
     208
     209        /*Output*/
     210        IssmDouble Dbar;
     211        inputs->GetInputAverage(&Dbar,DamageDbarEnum);
     212        return Dbar;
     213}
     214/*}}}*/
    177215/*FUNCTION Matice::GetVectorFromInputs{{{*/
    178216void  Matice::GetVectorFromInputs(Vector<IssmDouble>* vector,int input_enum){
    … …  
    210248void  Matice::GetViscosity2d(IssmDouble* pviscosity, IssmDouble* epsilon){
    211249        /*From a string tensor and a material object, return viscosity, using Glen's flow law.
    212                                                                                                     B
     250                                                                                                    (1-D) B
    213251          viscosity= -------------------------------------------------------------------
    214252                                                  2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
    … …  
    229267        /*Intermediary: */
    230268        IssmDouble A,e;
    231         IssmDouble B,n;
     269        IssmDouble B,D,n;
    232270
    233271        /*Get B and n*/
    234272        B=GetBbar();
    235273        n=GetN();
     274        D=GetDbar();
    236275
    237276        if (n==1){
    238277                /*Viscous behaviour! viscosity=B: */
    239                 viscosity=B/2;
     278                viscosity=(1-D)*B/2;
    240279        }
    241280        else{
    … …  
    257296                        else{
    258297                                e=(n-1)/(2*n);
    259                                 viscosity=B/(2*pow(A,e));
     298                                viscosity=(1-D)*B/(2*pow(A,e));
    260299                        }
    261300                }
    … …  
    266305        _assert_(B>0);
    267306        _assert_(n>0);
     307        _assert_(D>=0 && D<1);
    268308
    269309        /*Return: */
    … …  
    276316        /*Return viscosity accounting for steady state power law creep [Thomas and SSA, 1982]:
    277317         *
    278          *                                               B
     318         *                                               (1-D)*B
    279319         * viscosity3d= -------------------------------------------------------------------
    280320         *                      2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
    … …  
    295335        /*Intermediaries: */
    296336        IssmDouble A,e;
    297         IssmDouble B,n;
     337        IssmDouble B,D,n;
    298338
    299339        /*Get B and n*/
    300340        B=GetB();
     341        D=GetD();
    301342        n=GetN();
    302343
    303344        if (n==1){
    304345                /*Viscous behaviour! viscosity3d=B: */
    305                 viscosity3d=B/2;
     346                viscosity3d=(1-D)*B/2;
    306347        }
    307348        else{
    … …  
    328369                                e=(n-1)/2/n;
    329370
    330                                 viscosity3d=B/(2*pow(A,e));
     371                                viscosity3d=(1-D)*B/(2*pow(A,e));
    331372                        }
    332373                }
    … …  
    337378        _assert_(B>0);
    338379        _assert_(n>0);
     380        _assert_(D>=0 && D<1);
    339381
    340382        /*Assign output pointers:*/
    … …  
    346388        /*Return viscosity accounting for steady state power law creep [Thomas and SSA, 1982]:
    347389         *
    348          *                                          B
     390         *                                          (1-D)*B
    349391         * viscosity3d= -------------------------------------------------------------------
    350392         *                   2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
    … …  
    365407        /*Intermediaries: */
    366408        IssmDouble A,e;
    367         IssmDouble B,n;
     409        IssmDouble B,D,n;
    368410        IssmDouble eps0;
    369411
    … …  
    371413        eps0=pow(10.,-27);
    372414        B=GetB();
     415        D=GetD();
    373416        n=GetN();
    374417
    375418        if (n==1){
    376419                /*Viscous behaviour! viscosity3d=B: */
    377                 viscosity3d=B/2;
     420                viscosity3d=(1-D)*B/2;
    378421        }
    379422        else{
    … …  
    400443                        else{
    401444                                e=(n-1)/2/n;
    402                                 viscosity3d=B/(2*pow(A,e));
     445                                viscosity3d=(1-D)*B/(2*pow(A,e));
    403446                        }
    404447                }
    … …  
    409452        _assert_(B>0);
    410453        _assert_(n>0);
     454        _assert_(D>=0 && D<1);
    411455
    412456        /*Assign output pointers:*/
    … …  
    418462        /*Return viscosity accounting for steady state power law creep [Thomas and SSA, 1982]:
    419463         *
    420          *                                                                                              1
     464         *                                                                                              (1-D)
    421465         * viscosity= -------------------------------------------------------------------
    422466         *                                2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
    … …  
    434478        /*Intermediary value A and exponent e: */
    435479        IssmDouble A,e;
    436         IssmDouble B,n;
    437 
    438         /*Get B and n*/
    439         B=GetBbar();
     480        IssmDouble D,n;
     481
     482        /*Get D and n*/
     483        D=GetDbar();
    440484        n=GetN();
    441485
    … …  
    445489                exy=*(epsilon+2);
    446490
     491                /*Build viscosity: mu2=(1-D)/(2*A^e) */
     492                A=pow(exx,2)+pow(eyy,2)+pow(exy,2)+exx*eyy;
     493                if(A==0){
     494                        /*Maximum viscosity_complement for 0 shear areas: */
     495                        viscosity_complement=2.25*pow(10.,17);
     496                }
     497                else{
     498                        e=(n-1)/(2*n);
     499
     500                        viscosity_complement=(1-D)/(2*pow(A,e));
     501                }
     502        }
     503        else{
     504                viscosity_complement=4.5*pow(10.,17);
     505        }
     506
     507        /*Checks in debugging mode*/
     508        _assert_(D>=0 && D<1);
     509        _assert_(n>0);
     510        _assert_(viscosity_complement>0);
     511
     512        /*Return: */
     513        *pviscosity_complement=viscosity_complement;
     514}
     515/*}}}*/
     516/*FUNCTION Matice::GetViscosityDComplement {{{*/
     517void  Matice::GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* epsilon){
     518        /*Return viscosity derivative for control method d(mu)/dD:
     519         *
     520         *                                                                                             B
     521         * dviscosity= - -------------------------------------------------------------------
     522         *                                2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
     523         *
     524         * If epsilon is NULL, it means this is the first time Gradjb is being run, and we
     525         * return mu20, initial viscosity.
     526         */
     527
     528        /*output: */
     529        IssmDouble viscosity_complement;
     530
     531        /*input strain rate: */
     532        IssmDouble exx,eyy,exy;
     533
     534        /*Intermediary value A and exponent e: */
     535        IssmDouble A,e;
     536        IssmDouble B,n;
     537
     538        /*Get B and n*/
     539        B=GetBbar();
     540        n=GetN();
     541
     542        if(epsilon){
     543                exx=*(epsilon+0);
     544                eyy=*(epsilon+1);
     545                exy=*(epsilon+2);
     546
    447547                /*Build viscosity: mu2=B/(2*A^e) */
    448548                A=pow(exx,2)+pow(eyy,2)+pow(exy,2)+exx*eyy;
    449549                if(A==0){
    450550                        /*Maximum viscosity_complement for 0 shear areas: */
    451                         viscosity_complement=2.25*pow(10.,17);
     551                        viscosity_complement=- 2.25*pow(10.,17);
    452552                }
    453553                else{
    454554                        e=(n-1)/(2*n);
    455555
    456                         viscosity_complement=1/(2*pow(A,e));
     556                        viscosity_complement=- B/(2*pow(A,e));
    457557                }
    458558        }
    459559        else{
    460                 viscosity_complement=4.5*pow(10.,17);
     560                viscosity_complement=- 4.5*pow(10.,17);
    461561        }
    462562
    … …  
    464564        _assert_(B>0);
    465565        _assert_(n>0);
    466         _assert_(viscosity_complement>0);
     566        _assert_(viscosity_complement<0);
    467567
    468568        /*Return: */
    … …  
    729829                }
    730830
     831                /*Get D:*/
     832                if (iomodel->Data(DamageDEnum)) {
     833                        for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(DamageDEnum)[iomodel->elements[num_vertices*index+i]-1];
     834                        this->inputs->AddInput(new TriaInput(DamageDbarEnum,nodeinputs,P1Enum));
     835                }
     836
    731837                /*Control Inputs*/
    732838                #ifdef _HAVE_CONTROL_
    … …  
    743849                                                }
    744850                                                break;
     851                                        case DamageDbarEnum:
     852                                                if (iomodel->Data(DamageDEnum)){
     853                                                        _assert_(iomodel->Data(DamageDEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));
     854                                                        for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(DamageDEnum)[iomodel->elements[num_vertices*index+j]-1];
     855                                                        for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->Data(InversionMinParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
     856                                                        for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->Data(InversionMaxParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
     857                                                        this->inputs->AddInput(new ControlInput(DamageDbarEnum,TriaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
     858                                                }
     859                                                break;
    745860                                }
    746861                        }
    … …  
    770885                        this->inputs->AddInput(new PentaInput(MaterialsRheologyNEnum,nodeinputs,P1Enum));
    771886                }
    772 
     887               
     888                /*Get D*/
     889                if (iomodel->Data(DamageDEnum)) {
     890                        for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(DamageDEnum)[iomodel->elements[num_vertices*index+i]-1];
     891                        this->inputs->AddInput(new PentaInput(DamageDEnum,nodeinputs,P1Enum));
     892                }
     893               
    773894                /*Control Inputs*/
    774895                #ifdef _HAVE_CONTROL_
    … …  
    785906                                                }
    786907                                                break;
     908                                        case DamageDbarEnum:
     909                                                if (iomodel->Data(DamageDEnum)){
     910                                                        _assert_(iomodel->Data(DamageDEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));
     911                                                        for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(DamageDEnum)[iomodel->elements[num_vertices*index+j]-1];
     912                                                        for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->Data(InversionMinParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
     913                                                        for(j=0;j<num_vertices;j++)cmmaxinputs[j]=iomodel->Data(InversionMaxParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i];
     914                                                        this->inputs->AddInput(new ControlInput(DamageDEnum,PentaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));
     915                                                }
     916                                                break;
     917
    787918                                }
    788919                        }
    … …  
    804935                                name==MaterialsRheologyBbarEnum ||
    805936                                name==MaterialsRheologyNEnum ||
    806                                 name==QmuMaterialsRheologyBEnum
     937                                name==QmuMaterialsRheologyBEnum ||
     938                                name==DamageDEnum ||
     939                                name==DamageDbarEnum
     940
    807941                ){
    808942                return true;

  • issm/trunk-jpl/src/c/classes/Materials/Matice.h

    r16146 r16167  
    5757                void       GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon);
    5858                void       GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);
    59                 void       GetViscosityDComplement(IssmDouble*, IssmDouble*){_error_("not supported");};
     59                void       GetViscosityDComplement(IssmDouble*, IssmDouble*);
    6060                void       GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon);
    6161                void       GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon);
    … …  
    6464                IssmDouble GetB();
    6565                IssmDouble GetBbar();
    66                 IssmDouble GetZ(){_error_("not supported");};
    67                 IssmDouble GetD(){_error_("not supported");};
    68                 IssmDouble GetZbar(){_error_("not supported");};
    69                 IssmDouble GetDbar(){_error_("not supported");};
     66                IssmDouble GetZ();
     67                IssmDouble GetD();
     68                IssmDouble GetZbar();
     69                IssmDouble GetDbar();
    7070                IssmDouble GetN();
    7171                bool       IsInput(int name);

  • issm/trunk-jpl/src/c/classes/classes.h

    r15650 r16167  
    7777#include "./Materials/Material.h"
    7878#include "./Materials/Matice.h"
    79 #include "./Materials/Matdamageice.h"
    8079#include "./Materials/Matpar.h"
    8180

  • issm/trunk-jpl/src/c/modules/Gradjx/Gradjx.cpp

    r14999 r16167  
    3434
    3535                for(j=0;j<elements->Size();j++){
     36
    3637                        Element* element=(Element*)elements->GetObjectByOffset(j);
    3738                        element->Gradj(gradient_list[i],control_type[i],i);

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

    r16162 r16167  
    5454        switch(materials_type){
    5555                case MaticeEnum:
    56                         iomodel->FetchData(2,MaterialsRheologyBEnum,MaterialsRheologyNEnum);
     56                        iomodel->FetchData(3,MaterialsRheologyBEnum,MaterialsRheologyNEnum,DamageDEnum);
    5757                        for (i=0;i<iomodel->numberofelements;i++) if(iomodel->my_elements[i]) materials->AddObject(new Matice(i+1,i,iomodel));
    5858                        if(dakota_analysis){
    … …  
    6262           #endif
    6363                        }
    64                         break;
    65                 case MatdamageiceEnum:
    66                         iomodel->FetchData(3,MaterialsRheologyBEnum,MaterialsRheologyNEnum,DamageDEnum);
    67                         for (i=0;i<iomodel->numberofelements;i++) if(iomodel->my_elements[i]) materials->AddObject(new Matdamageice(i+1,i,iomodel));
    6864                        break;
    6965                default:

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Stressbalance/UpdateElementsStressbalance.cpp

    r16162 r16167  
    101101        iomodel->FetchDataToInput(elements,LoadingforceXEnum);
    102102        iomodel->FetchDataToInput(elements,LoadingforceYEnum);
    103         if(materials_type==MatdamageiceEnum){
    104                 iomodel->FetchDataToInput(elements,DamageDEnum);
    105         }
     103        iomodel->FetchDataToInput(elements,DamageDEnum);
     104
    106105        if(iomodel->dim==3){
    107106                iomodel->FetchDataToInput(elements,MeshElementonbedEnum);

  • issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

    r16162 r16167  
    152152        MaterialsRheologyLawEnum,
    153153        MaterialsRheologyNEnum,
    154         MaterialsRheologyZEnum,
    155         MaterialsRheologyZbarEnum,
    156154        DamageDEnum,
    157155        DamageDbarEnum,
    … …  
    373371        SSA3dIceFrontEnum,
    374372        MaticeEnum,
    375         MatdamageiceEnum,
    376373        MatparEnum,
    377374        NodeEnum,

  • issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

    r16162 r16167  
    160160                case MaterialsRheologyLawEnum : return "MaterialsRheologyLaw";
    161161                case MaterialsRheologyNEnum : return "MaterialsRheologyN";
    162                 case MaterialsRheologyZEnum : return "MaterialsRheologyZ";
    163                 case MaterialsRheologyZbarEnum : return "MaterialsRheologyZbar";
    164162                case DamageDEnum : return "DamageD";
    165163                case DamageDbarEnum : return "DamageDbar";
    … …  
    373371                case SSA3dIceFrontEnum : return "SSA3dIceFront";
    374372                case MaticeEnum : return "Matice";
    375                 case MatdamageiceEnum : return "Matdamageice";
    376373                case MatparEnum : return "Matpar";
    377374                case NodeEnum : return "Node";

  • issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

    r16162 r16167  
    163163              else if (strcmp(name,"MaterialsRheologyLaw")==0) return MaterialsRheologyLawEnum;
    164164              else if (strcmp(name,"MaterialsRheologyN")==0) return MaterialsRheologyNEnum;
    165               else if (strcmp(name,"MaterialsRheologyZ")==0) return MaterialsRheologyZEnum;
    166               else if (strcmp(name,"MaterialsRheologyZbar")==0) return MaterialsRheologyZbarEnum;
    167165              else if (strcmp(name,"DamageD")==0) return DamageDEnum;
    168166              else if (strcmp(name,"DamageDbar")==0) return DamageDbarEnum;
    … …  
    260258              else if (strcmp(name,"ThermalMaxiter")==0) return ThermalMaxiterEnum;
    261259              else if (strcmp(name,"ThermalPenaltyFactor")==0) return ThermalPenaltyFactorEnum;
     260              else if (strcmp(name,"ThermalPenaltyLock")==0) return ThermalPenaltyLockEnum;
     261              else if (strcmp(name,"ThermalPenaltyThreshold")==0) return ThermalPenaltyThresholdEnum;
    262262         else stage=3;
    263263   }
    264264   if(stage==3){
    265               if (strcmp(name,"ThermalPenaltyLock")==0) return ThermalPenaltyLockEnum;
    266               else if (strcmp(name,"ThermalPenaltyThreshold")==0) return ThermalPenaltyThresholdEnum;
    267               else if (strcmp(name,"ThermalSpctemperature")==0) return ThermalSpctemperatureEnum;
     265              if (strcmp(name,"ThermalSpctemperature")==0) return ThermalSpctemperatureEnum;
    268266              else if (strcmp(name,"ThermalStabilization")==0) return ThermalStabilizationEnum;
    269267              else if (strcmp(name,"GiaMantleViscosity")==0) return GiaMantleViscosityEnum;
    … …  
    379377              else if (strcmp(name,"SSA3dIceFront")==0) return SSA3dIceFrontEnum;
    380378              else if (strcmp(name,"Matice")==0) return MaticeEnum;
    381               else if (strcmp(name,"Matdamageice")==0) return MatdamageiceEnum;
    382379              else if (strcmp(name,"Matpar")==0) return MatparEnum;
    383380              else if (strcmp(name,"Node")==0) return NodeEnum;
    384381              else if (strcmp(name,"Numericalflux")==0) return NumericalfluxEnum;
     382              else if (strcmp(name,"NumericalfluxType")==0) return NumericalfluxTypeEnum;
     383              else if (strcmp(name,"Param")==0) return ParamEnum;
     384              else if (strcmp(name,"L1L2IceFront")==0) return L1L2IceFrontEnum;
    385385         else stage=4;
    386386   }
    387387   if(stage==4){
    388               if (strcmp(name,"NumericalfluxType")==0) return NumericalfluxTypeEnum;
    389               else if (strcmp(name,"Param")==0) return ParamEnum;
    390               else if (strcmp(name,"L1L2IceFront")==0) return L1L2IceFrontEnum;
    391               else if (strcmp(name,"HOIceFront")==0) return HOIceFrontEnum;
     388              if (strcmp(name,"HOIceFront")==0) return HOIceFrontEnum;
    392389              else if (strcmp(name,"Pengrid")==0) return PengridEnum;
    393390              else if (strcmp(name,"Penpair")==0) return PenpairEnum;
    … …  
    506503              else if (strcmp(name,"StressTensorxz")==0) return StressTensorxzEnum;
    507504              else if (strcmp(name,"StressTensoryy")==0) return StressTensoryyEnum;
     505              else if (strcmp(name,"StressTensoryz")==0) return StressTensoryzEnum;
     506              else if (strcmp(name,"StressTensorzz")==0) return StressTensorzzEnum;
     507              else if (strcmp(name,"GiaCrossSectionShape")==0) return GiaCrossSectionShapeEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"StressTensoryz")==0) return StressTensoryzEnum;
    512               else if (strcmp(name,"StressTensorzz")==0) return StressTensorzzEnum;
    513               else if (strcmp(name,"GiaCrossSectionShape")==0) return GiaCrossSectionShapeEnum;
    514               else if (strcmp(name,"GiadWdt")==0) return GiadWdtEnum;
     511              if (strcmp(name,"GiadWdt")==0) return GiadWdtEnum;
    515512              else if (strcmp(name,"GiaW")==0) return GiaWEnum;
    516513              else if (strcmp(name,"P0")==0) return P0Enum;

  • issm/trunk-jpl/src/m/classes/damage.m

    r16162 r16167  
    5656                function md = checkconsistency(obj,md,solution,analyses) % {{{
    5757                       
     58                        md = checkfield(md,'damage.D','>=0',0,'size',[md.mesh.numberofvertices 1]);
    5859                        if obj.isdamage,
    59                                 md = checkfield(md,'damage.D','>=0',0,'size',[md.mesh.numberofvertices 1]);
    6060                                md = checkfield(md,'damage.law','values',{'undamaged','pralong'});
    6161                        end
    … …  
    8686               
    8787                        WriteData(fid,'object',obj,'class','damage','fieldname','isdamage','format','DoubleMat','mattype',1);
     88                        WriteData(fid,'object',obj,'class','damage','fieldname','D','format','DoubleMat','mattype',1);
    8889                        if obj.isdamage,
    89                                 WriteData(fid,'object',obj,'class','damage','fieldname','D','format','DoubleMat','mattype',1);
    9090                                WriteData(fid,'object',obj,'class','damage','fieldname','law','format','String');
    9191                        end

  • issm/trunk-jpl/src/m/classes/inversion.m

    r16148 r16167  
    116116                        md = checkfield(md,'inversion.incomplete_adjoint','values',[0 1]);
    117117                        md = checkfield(md,'inversion.control_parameters','cell',1,'values',...
    118                                 {'BalancethicknessThickeningRate' 'FrictionCoefficient' 'MaterialsRheologyBbar' 'MaterialsDamageDbar' 'Vx' 'Vy' 'Thickness'});
     118                                {'BalancethicknessThickeningRate' 'FrictionCoefficient' 'MaterialsRheologyBbar' 'DamageDbar' '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

    r16148 r16167  
    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','MaterialsDamageDbar','Vx','Vy'])
     126                md = checkfield(md,'inversion.control_parameters','cell',1,'values',['BalancethicknessThickeningRate','FrictionCoefficient','MaterialsRheologyBbar','DamageDbar','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/enum/EnumDefinitions.py

    r16162 r16167  
    152152def MaterialsRheologyLawEnum(): return StringToEnum("MaterialsRheologyLaw")[0]
    153153def MaterialsRheologyNEnum(): return StringToEnum("MaterialsRheologyN")[0]
    154 def MaterialsRheologyZEnum(): return StringToEnum("MaterialsRheologyZ")[0]
    155 def MaterialsRheologyZbarEnum(): return StringToEnum("MaterialsRheologyZbar")[0]
    156154def DamageDEnum(): return StringToEnum("DamageD")[0]
    157155def DamageDbarEnum(): return StringToEnum("DamageDbar")[0]
    … …  
    365363def SSA3dIceFrontEnum(): return StringToEnum("SSA3dIceFront")[0]
    366364def MaticeEnum(): return StringToEnum("Matice")[0]
    367 def MatdamageiceEnum(): return StringToEnum("Matdamageice")[0]
    368365def MatparEnum(): return StringToEnum("Matpar")[0]
    369366def NodeEnum(): return StringToEnum("Node")[0]

  • issm/trunk-jpl/test/NightlyRun/test213.m

    r15771 r16167  
    1818md.inversion.vx_obs=md.initialization.vx; md.inversion.vy_obs=md.initialization.vy;
    1919
    20 md.cluster=generic('name',oshostname(),'np',3);
     20md.verbose.control=true;
     21md.verbose.solution=true;
     22
     23md.cluster=generic('name',oshostname(),'np',1);
    2124md=solve(md,StressbalanceSolutionEnum());
    2225

  • issm/trunk-jpl/test/NightlyRun/test270.m

    r16161 r16167  
    22md=setmask(md,'all','');
    33md=parameterize(md,'../Par/SquareShelf.par');
    4 md.materials=matdamageice();
    5 md.materials.rheology_B=paterson(md.initialization.temperature);
    6 md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    7 md.damage.isdamage=1;
    84md.damage.D=0.5*ones(md.mesh.numberofvertices,1);
    95md=setflowequation(md,'SSA','all');

  • issm/trunk-jpl/test/NightlyRun/test270.py

    r16154 r16167  
    44from setmask import *
    55from parameterize import *
    6 from matdamageice import *
    76from paterson import *
    87from setflowequation import *
    … …  
    1413md=setmask(md,'all','')
    1514md=parameterize(md,'../Par/SquareShelf.py')
    16 md.materials=matdamageice()
    17 md.materials.rheology_B=paterson(md.initialization.temperature)
    18 md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    19 md.materials.damage_D=0.5*numpy.ones((md.mesh.numberofvertices,1))
     15md.damage.D=0.5*numpy.ones((md.mesh.numberofvertices,1))
    2016md=setflowequation(md,'SSA','all')
    2117md.cluster=generic('name',oshostname(),'np',3)

  • issm/trunk-jpl/test/NightlyRun/test272.m

    r16161 r16167  
    22md=setmask(md,'all','');
    33md=parameterize(md,'../Par/SquareShelf.par');
    4 md.materials=matdamageice();
    5 md.materials.rheology_B=paterson(md.initialization.temperature);
    6 md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    7 md.damage.isdamage=1;
    84md.damage.D=0.5*ones(md.mesh.numberofvertices,1);
    95md=setflowequation(md,'SSA','all');
    … …  
    128%control parameters
    139md.inversion.iscontrol=1;
    14 md.inversion.control_parameters={'MaterialsDamageDbar'};
     10md.inversion.control_parameters={'DamageDbar'};
    1511md.inversion.min_parameters=zeros(md.mesh.numberofvertices,1);
    1612md.inversion.max_parameters=(1-10^-13)*ones(md.mesh.numberofvertices,1);
    … …  
    2824
    2925%Fields and tolerances to track changes
    30 field_names     ={'Gradient','Misfits','MaterialsDamageDbar','Pressure','Vel','Vx','Vy'};
     26field_names     ={'Gradient','Misfits','DamageDbar','Pressure','Vel','Vx','Vy'};
    3127field_tolerances={1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12};
    3228field_values={...
    3329   (md.results.StressbalanceSolution.Gradient1),...
    3430   (md.results.StressbalanceSolution.J),...
    35    (md.results.StressbalanceSolution.MaterialsDamageDbar),...
     31   (md.results.StressbalanceSolution.DamageDbar),...
    3632   (md.results.StressbalanceSolution.Pressure),...
    3733   (md.results.StressbalanceSolution.Vel),...

  • issm/trunk-jpl/test/NightlyRun/test272.py

    r16154 r16167  
    44from setmask import *
    55from parameterize import *
    6 from matdamageice import *
    76from paterson import *
    87from setflowequation import *
    … …  
    1413md=setmask(md,'all','')
    1514md=parameterize(md,'../Par/SquareShelf.py')
    16 md.materials=matdamageice()
    17 md.materials.rheology_B=paterson(md.initialization.temperature)
    18 md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    19 md.materials.damage_D=0.5*numpy.ones((md.mesh.numberofvertices,1))
     15md.damage.D=0.5*numpy.ones((md.mesh.numberofvertices,1))
    2016md=setflowequation(md,'SSA','all')
    2117
    2218#control parameters
    2319md.inversion.iscontrol=1
    24 md.inversion.control_parameters=['MaterialsDamageDbar']
     20md.inversion.control_parameters=['DamageDbar']
    2521md.inversion.min_parameters=10**-13*numpy.ones((md.mesh.numberofvertices,1))
    2622md.inversion.max_parameters=numpy.ones((md.mesh.numberofvertices,1))
    … …  
    3834
    3935#Fields and tolerances to track changes
    40 field_names     =['Gradient','Misfits','MaterialsRheologyDbar','Pressure','Vel','Vx','Vy']
     36field_names     =['Gradient','Misfits','DamageDbar','Pressure','Vel','Vx','Vy']
    4137field_tolerances=[1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12]
    4238field_values=[\
    4339   md.results.StressbalanceSolution.Gradient1,\
    4440   md.results.StressbalanceSolution.J,\
    45    md.results.StressbalanceSolution.MaterialsDamageDbar,\
     41   md.results.StressbalanceSolution.DamageDbar,\
    4642   md.results.StressbalanceSolution.Pressure,\
    4743   md.results.StressbalanceSolution.Vel,\

  • issm/trunk-jpl/test/NightlyRun/test274.m

    r16161 r16167  
    33md=setmask(md,'all','');
    44md=parameterize(md,'../Par/SquareShelf2.par');
    5 md.materials=matdamageice();
    6 md.materials.rheology_B=paterson(md.initialization.temperature);
    7 md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    8 md.damage.isdamage=1;
    95md.damage.D=0.5*ones(md.mesh.numberofvertices,1);
    106md=setflowequation(md,'SSA','all');

  • issm/trunk-jpl/test/NightlyRun/test274.py

    r16154 r16167  
    55from setmask import *
    66from parameterize import *
    7 from matdamageice import *
    87from paterson import *
    98from setflowequation import *
    … …  
    1615md=setmask(md,'all','')
    1716md=parameterize(md,'../Par/SquareShelf2.py')
    18 md.materials=matdamageice()
    19 md.materials.rheology_B=paterson(md.initialization.temperature)
    20 md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    21 md.materials.damage_D=0.5*numpy.ones((md.mesh.numberofvertices,1))
     17md.damage.D=0.5*numpy.ones((md.mesh.numberofvertices,1))
    2218md=setflowequation(md,'SSA','all')
    2319

  • issm/trunk-jpl/test/Par/79North.par

    r15988 r16167  
    2020md.materials.rheology_B=paterson(md.initialization.temperature);
    2121md.materials.rheology_n=3*ones(md.mesh.numberofelements,1);
     22
     23%Damage
     24md.damage.D=zeros(md.mesh.numberofvertices,1);
    2225
    2326%Friction

  • issm/trunk-jpl/test/Par/79North.py

    r15988 r16167  
    3333md.initialization.temperature=md.initialization.temperature
    3434
     35#Damage
     36md.damage.D=zeros((md.mesh.numberofvertices,1))
     37
    3538#Friction
    3639md.friction.coefficient=50.*ones((md.mesh.numberofvertices,1))

  • issm/trunk-jpl/test/Par/GiaBenchmarksAB.par

    r15988 r16167  
    4444md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    4545
     46%Damage
     47md.damage.D=zeros(md.mesh.numberofvertices,1);
     48
    4649%Friction
    4750md.friction.coefficient=20.*ones(md.mesh.numberofvertices,1);

  • issm/trunk-jpl/test/Par/GiaBenchmarksCD.par

    r15988 r16167  
    4343md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    4444
     45%Damage
     46md.damage.D=zeros(md.mesh.numberofvertices,1);
     47
    4548%Friction
    4649md.friction.coefficient=20.*ones(md.mesh.numberofvertices,1);

  • issm/trunk-jpl/test/Par/ISMIPA.par

    r15988 r16167  
    1616md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    1717
     18%Damage
     19md.damage.D=zeros(md.mesh.numberofvertices,1);
     20
    1821disp('      boundary conditions for stressbalance model');
    1922%Create node on boundary first (because we cannot use mesh)

  • issm/trunk-jpl/test/Par/ISMIPA.py

    r15988 r16167  
    1919md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    2020
     21#Damage
     22md.damage.D=zeros((md.mesh.numberofvertices,1))
     23
    2124print "      boundary conditions for stressbalance model"
    2225#Create node on boundary first (because we cannot use mesh)

  • issm/trunk-jpl/test/Par/ISMIPB.par

    r15988 r16167  
    1616md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    1717
     18%Damage
     19md.damage.D=zeros(md.mesh.numberofvertices,1);
     20
    1821disp('      boundary conditions for stressbalance model');
    1922%Create node on boundary first (because we cannot use mesh)

  • issm/trunk-jpl/test/Par/ISMIPB.py

    r15988 r16167  
    1919md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    2020
     21#Damage
     22md.damage.D=zeros((md.mesh.numberofvertices,1))
     23
    2124print "      boundary conditions for stressbalance model"
    2225#Create node on boundary first (because we cannot use mesh)

  • issm/trunk-jpl/test/Par/ISMIPC.par

    r15988 r16167  
    1717md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    1818
     19%Damage
     20md.damage.D=zeros(md.mesh.numberofvertices,1);
     21
    1922disp('      boundary conditions for stressbalance model:');
    2023%Create node on boundary first (because we can not use mesh)

  • issm/trunk-jpl/test/Par/ISMIPC.py

    r15988 r16167  
    2020md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    2121
     22#Damage
     23md.damage.D=zeros((md.mesh.numberofvertices,1))
     24
    2225print "      boundary conditions for stressbalance model:"
    2326#Create node on boundary first (because we can not use mesh)

  • issm/trunk-jpl/test/Par/ISMIPD.par

    r15988 r16167  
    1616md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    1717
     18%Damage
     19md.damage.D=zeros(md.mesh.numberofvertices,1);
     20
    1821disp('      boundary conditions for stressbalance model:');
    1922%Create node on boundary first (because we can not use mesh)

  • issm/trunk-jpl/test/Par/ISMIPD.py

    r15988 r16167  
    1919md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    2020
     21#Damage
     22md.damage.D=zeros((md.mesh.numberofvertices,1))
     23
    2124print "      boundary conditions for stressbalance model:"
    2225#Create node on boundary first (because we can not use mesh)

  • issm/trunk-jpl/test/Par/ISMIPE.par

    r15771 r16167  
    2626md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    2727
     28%Damage
     29md.damage.D=zeros(md.mesh.numberofvertices,1);
     30
    2831disp('      boundary conditions for stressbalance model:');
    2932%Create node on boundary first (because we can not use mesh)

  • issm/trunk-jpl/test/Par/ISMIPE.py

    r15771 r16167  
    3131md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    3232
     33#Damage
     34md.damage.D=zeros((md.mesh.numberofvertices,1))
     35
    3336print "      boundary conditions for stressbalance model:"
    3437#Create node on boundary first (because we can not use mesh)

  • issm/trunk-jpl/test/Par/ISMIPF.par

    r15771 r16167  
    1717md.materials.rheology_n=1.*ones(md.mesh.numberofelements,1);
    1818md.materials.rheology_law='None';
     19
     20%Damage
     21md.damage.D=zeros(md.mesh.numberofvertices,1);
    1922
    2023disp('      boundary conditions for stressbalance model');

  • issm/trunk-jpl/test/Par/ISMIPF.py

    r15771 r16167  
    2020md.materials.rheology_n=1.*numpy.ones((md.mesh.numberofelements,1))
    2121md.materials.rheology_law='None'
     22
     23#Damage
     24md.damage.D=zeros((md.mesh.numberofvertices,1))
    2225
    2326print "      boundary conditions for stressbalance model"

  • issm/trunk-jpl/test/Par/Pig.par

    r15988 r16167  
    2626md.initialization.temperature=md.initialization.temperature;
    2727
     28%Damage
     29md.damage.D=zeros(md.mesh.numberofvertices,1);
     30
    2831%Friction
    2932md.friction.coefficient=50*ones(md.mesh.numberofvertices,1);

  • issm/trunk-jpl/test/Par/Pig.py

    r15988 r16167  
    3737md.initialization.temperature=md.initialization.temperature
    3838
     39#Damage
     40md.damage.D=zeros((md.mesh.numberofvertices,1))
     41
    3942#Friction
    4043md.friction.coefficient=50.*ones((md.mesh.numberofvertices,1))

  • issm/trunk-jpl/test/Par/RoundSheetEISMINT.par

    r15988 r16167  
    2020md.materials.rheology_B=6.81*10^7*ones(md.mesh.numberofvertices,1); %to have the same B as the analytical solution
    2121md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
     22
     23%Damage
     24md.damage.D=zeros(md.mesh.numberofvertices,1);
    2225
    2326disp('      creating surface mass balance');

  • issm/trunk-jpl/test/Par/RoundSheetEISMINT.py

    r15988 r16167  
    2323md.materials.rheology_B=6.81*10**7*numpy.ones((md.mesh.numberofvertices,1))    #to have the same B as the analytical solution
    2424md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
     25
     26#Damage
     27md.damage.D=zeros((md.mesh.numberofvertices,1))
    2528
    2629print "      creating surface mass balance"

  • issm/trunk-jpl/test/Par/RoundSheetShelf.par

    r15994 r16167  
    6161md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    6262
     63%Damage
     64md.damage.D=zeros(md.mesh.numberofvertices,1);
     65
    6366%Surface mass balance and basal melting
    6467md.surfaceforcings.mass_balance=-10.*ones(md.mesh.numberofvertices,1);

  • issm/trunk-jpl/test/Par/RoundSheetShelf.py

    r16005 r16167  
    6868md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    6969
     70#Damage
     71md.damage.D=zeros((md.mesh.numberofvertices,1))
     72
    7073#Surface mass balance and basal melting
    7174md.surfaceforcings.mass_balance=-10.*numpy.ones((md.mesh.numberofvertices,1))

  • issm/trunk-jpl/test/Par/RoundSheetStaticEISMINT.par

    r15988 r16167  
    2424md.materials.rheology_B=6.81*10^7*ones(md.mesh.numberofvertices,1); %to have the same B as the analytical solution
    2525md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
     26
     27%Damage
     28md.damage.D=zeros(md.mesh.numberofvertices,1);
    2629
    2730disp('      creating surface mass balance');

  • issm/trunk-jpl/test/Par/RoundSheetStaticEISMINT.py

    r15988 r16167  
    2828md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    2929
     30#Damage
     31md.damage.D=zeros((md.mesh.numberofvertices,1))
     32
    3033print "      creating surface mass balance"
    3134smb_max=0.5    #m/yr

  • issm/trunk-jpl/test/Par/SquareEISMINT.par

    r15988 r16167  
    2626md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    2727
     28%Damage
     29md.damage.D=zeros(md.mesh.numberofvertices,1);
     30
    2831disp('      creating surface mass balance');
    2932md.surfaceforcings.mass_balance=0.2*ones(md.mesh.numberofvertices,1); %0m/a

  • issm/trunk-jpl/test/Par/SquareEISMINT.py

    r15988 r16167  
    2929md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    3030
     31#Damage
     32md.damage.D=zeros((md.mesh.numberofvertices,1))
     33
    3134print "      creating surface mass balance"
    3235md.surfaceforcings.mass_balance=0.2*numpy.ones((md.mesh.numberofvertices,1))    #0m/a

  • issm/trunk-jpl/test/Par/SquareSheetConstrained.par

    r16063 r16167  
    2929md.materials.rheology_B=paterson(md.initialization.temperature);
    3030md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
     31
     32%Damage
     33md.damage.D=zeros(md.mesh.numberofvertices,1);
    3134
    3235%Friction

  • issm/trunk-jpl/test/Par/SquareSheetConstrained.py

    r16063 r16167  
    4040md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    4141
     42#Damage
     43md.damage.D=zeros((md.mesh.numberofvertices,1))
     44
    4245#Friction
    4346md.friction.coefficient=20.*numpy.ones((md.mesh.numberofvertices,1))

  • issm/trunk-jpl/test/Par/SquareSheetShelf.par

    r16063 r16167  
    3232md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    3333
     34%Damage
     35md.damage.D=zeros(md.mesh.numberofvertices,1);
     36
    3437%Accumulation and melting
    3538md.surfaceforcings.mass_balance=10.*ones(md.mesh.numberofvertices,1);

  • issm/trunk-jpl/test/Par/SquareSheetShelf.py

    r16063 r16167  
    4343md.materials.rheology_n=3.*ones((md.mesh.numberofelements,1))
    4444
     45#Damage
     46md.damage.D=zeros((md.mesh.numberofvertices,1))
     47
    4548#Accumulation and melting
    4649md.surfaceforcings.mass_balance=10.*numpy.ones((md.mesh.numberofvertices,1))

  • issm/trunk-jpl/test/Par/SquareShelf.par

    r16063 r16167  
    2929md.materials.rheology_n=3*ones(md.mesh.numberofelements,1);
    3030
     31%Damage
     32md.damage.D=zeros(md.mesh.numberofvertices,1);
     33
    3134%Friction
    3235md.friction.coefficient=20*ones(md.mesh.numberofvertices,1);

  • issm/trunk-jpl/test/Par/SquareShelf.py

    r16063 r16167  
    5959md.materials.rheology_n = 3.*ones((md.mesh.numberofelements,1))
    6060
     61#Damage
     62md.damage.D=zeros((md.mesh.numberofvertices,1))
     63
    6164#Friction
    6265md.friction.coefficient = 20.*ones((md.mesh.numberofvertices,1))

  • issm/trunk-jpl/test/Par/SquareShelf2.par

    r16063 r16167  
    2929md.materials.rheology_n=3*ones(md.mesh.numberofelements,1);
    3030
     31%Damage
     32md.damage.D=zeros(md.mesh.numberofvertices,1);
     33
    3134%Friction
    3235md.friction.coefficient=20*ones(md.mesh.numberofvertices,1);

  • issm/trunk-jpl/test/Par/SquareShelf2.py

    r16063 r16167  
    5959md.materials.rheology_n = 3.*ones((md.mesh.numberofelements,1))
    6060
     61#Damage
     62md.damage.D=zeros((md.mesh.numberofvertices,1))
     63
    6164#Friction
    6265md.friction.coefficient = 20.*ones((md.mesh.numberofvertices,1))

  • issm/trunk-jpl/test/Par/SquareShelfConstrained.par

    r16063 r16167  
    2929md.materials.rheology_n=3*ones(md.mesh.numberofelements,1);
    3030
     31%Damage
     32md.damage.D=zeros(md.mesh.numberofvertices,1);
     33
    3134%Surface mass balance and basal melting
    3235md.surfaceforcings.mass_balance=10*ones(md.mesh.numberofvertices,1);

  • issm/trunk-jpl/test/Par/SquareShelfConstrained.py

    r16063 r16167  
    4040md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    4141
     42#Damage
     43md.damage.D=zeros((md.mesh.numberofvertices,1))
     44
    4245#Surface mass balance and basal melting
    4346md.surfaceforcings.mass_balance=10.*numpy.ones((md.mesh.numberofvertices,1))

  • issm/trunk-jpl/test/Par/SquareThermal.par

    r15988 r16167  
    2727md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1);
    2828
     29%Damage
     30md.damage.D=zeros(md.mesh.numberofvertices,1);
     31
    2932disp('      creating surface mass balance');
    3033md.surfaceforcings.mass_balance=ones(md.mesh.numberofvertices,1)/md.constants.yts; %1m/a

  • issm/trunk-jpl/test/Par/SquareThermal.py

    r15988 r16167  
    3131md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1))
    3232
     33#Damage
     34md.damage.D=zeros((md.mesh.numberofvertices,1))
     35
    3336print "      creating surface mass balance"
    3437md.surfaceforcings.mass_balance=numpy.ones((md.mesh.numberofvertices,1))/md.constants.yts    #1m/a
Note: See TracChangeset for help on using the changeset viewer.