Changeset 16146
- Timestamp:
- 09/17/13 12:44:19 (12 years ago)
- Location:
- issm/trunk-jpl
- Files:
-
- 2 added
- 25 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk-jpl/src/c/classes/Elements/Tria.cpp
r16145 r16146 1006 1006 /*Computeportion of the element that is grounded*/ 1007 1007 1008 int normal_orientation =0;1008 int normal_orientation; 1009 1009 IssmDouble s1,s2; 1010 1010 IssmDouble levelset[NUMVERTICES]; … … 1013 1013 GetInputListOnVertices(&levelset[0],levelsetenum); 1014 1014 1015 if(levelset[0]*levelset[1]>0 .){ //Nodes 0 and 1 are similar, so points must be found on segment 0-2 and 1-21015 if(levelset[0]*levelset[1]>0){ //Nodes 0 and 1 are similar, so points must be found on segment 0-2 and 1-2 1016 1016 /*Portion of the segments*/ 1017 1017 s1=levelset[2]/(levelset[2]-levelset[1]); 1018 1018 s2=levelset[2]/(levelset[2]-levelset[0]); 1019 1019 1020 if(levelset[2]>0 .) normal_orientation=1;1020 if(levelset[2]>0) normal_orientation=1; 1021 1021 /*New point 1*/ 1022 1022 xyz_zero[3*normal_orientation+0]=xyz_list[2][0]+s1*(xyz_list[1][0]-xyz_list[2][0]); … … 1029 1029 xyz_zero[3*(1-normal_orientation)+2]=xyz_list[2][2]+s2*(xyz_list[0][2]-xyz_list[2][2]); 1030 1030 } 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-21031 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 1032 1032 /*Portion of the segments*/ 1033 1033 s1=levelset[0]/(levelset[0]-levelset[2]); 1034 1034 s2=levelset[0]/(levelset[0]-levelset[1]); 1035 1035 1036 if(levelset[0]>0 .) normal_orientation=1;1036 if(levelset[0]>0) normal_orientation=1; 1037 1037 /*New point 1*/ 1038 1038 xyz_zero[3*normal_orientation+0]=xyz_list[0][0]+s1*(xyz_list[2][0]-xyz_list[0][0]); … … 1045 1045 xyz_zero[3*(1-normal_orientation)+2]=xyz_list[0][2]+s2*(xyz_list[1][2]-xyz_list[0][2]); 1046 1046 } 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-21047 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 1048 1048 /*Portion of the segments*/ 1049 1049 s1=levelset[1]/(levelset[1]-levelset[0]); 1050 1050 s2=levelset[1]/(levelset[1]-levelset[2]); 1051 1051 1052 if(levelset[1]>0 .) normal_orientation=1;1052 if(levelset[1]>0) normal_orientation=1; 1053 1053 /*New point 0*/ 1054 1054 xyz_zero[3*normal_orientation+0]=xyz_list[1][0]+s1*(xyz_list[0][0]-xyz_list[1][0]); … … 1061 1061 xyz_zero[3*(1-normal_orientation)+2]=xyz_list[1][2]+s2*(xyz_list[2][2]-xyz_list[1][2]); 1062 1062 } 1063 else if(levelset[0]==0 . && levelset[1]==0.){ //front is on point 0 and 11063 else if(levelset[0]==0 && levelset[1]==0){ //front is on point 0 and 1 1064 1064 xyz_zero[3*0+0]=xyz_list[0][0]; 1065 1065 xyz_zero[3*0+1]=xyz_list[0][1]; … … 1071 1071 xyz_zero[3*1+2]=xyz_list[1][2]; 1072 1072 } 1073 else if(levelset[0]==0 . && levelset[2]==0.){ //front is on point 0 and 11073 else if(levelset[0]==0 && levelset[2]==0){ //front is on point 0 and 1 1074 1074 xyz_zero[3*0+0]=xyz_list[2][0]; 1075 1075 xyz_zero[3*0+1]=xyz_list[2][1]; … … 1081 1081 xyz_zero[3*1+2]=xyz_list[0][2]; 1082 1082 } 1083 else if(levelset[1]==0 . && levelset[2]==0.){ //front is on point 0 and 11083 else if(levelset[1]==0 && levelset[2]==0){ //front is on point 0 and 1 1084 1084 xyz_zero[3*0+0]=xyz_list[1][0]; 1085 1085 xyz_zero[3*0+1]=xyz_list[1][1]; … … 1451 1451 1452 1452 /*Go through all the input objects, and find the one corresponding to enum_type, if it exists: */ 1453 if (enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheology ZbarEnum)1453 if (enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyDbarEnum) 1454 1454 input=this->material->inputs->GetInput(enum_type); 1455 1455 else … … 1567 1567 break; 1568 1568 case MaterialsRheologyBbarEnum: 1569 case MaterialsRheology ZbarEnum:1569 case MaterialsRheologyDbarEnum: 1570 1570 /*Material will take care of it*/ break; 1571 1571 default: … … 1767 1767 void Tria::InputUpdateFromVector(IssmDouble* vector, int name, int type){ 1768 1768 1769 1769 1770 /*Check that name is an element input*/ 1770 1771 if (!IsInput(name)) return; … … 1780 1781 } 1781 1782 /*update input*/ 1782 if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheology ZEnum || name==MaterialsRheologyZbarEnum){1783 if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheologyDEnum || name==MaterialsRheologyDbarEnum){ 1783 1784 material->inputs->AddInput(new TriaInput(name,values,P1Enum)); 1784 1785 } … … 1793 1794 } 1794 1795 /*update input*/ 1795 if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheology ZEnum || name==MaterialsRheologyZbarEnum){1796 if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsRheologyDEnum || name==MaterialsRheologyDbarEnum){ 1796 1797 material->inputs->AddInput(new TriaInput(name,values,P1Enum)); 1797 1798 } … … 1961 1962 name==FrictionCoefficientEnum || 1962 1963 name==MaterialsRheologyBbarEnum || 1963 name==MaterialsRheology ZbarEnum ||1964 name==MaterialsRheologyDbarEnum || 1964 1965 name==GradientEnum || 1965 1966 name==OldGradientEnum || … … 2793 2794 *presponse=this->material->GetBbar(); 2794 2795 break; 2795 case MaterialsRheologyZbarEnum: 2796 *presponse=this->material->GetZbar(); 2797 break; 2796 2798 2797 case VelEnum:{ 2799 2798 … … 3731 3730 for(int i=0;i<num_controls;i++){ 3732 3731 3733 if(control_type[i]==MaterialsRheologyBbarEnum || control_type[i]==MaterialsRheology ZbarEnum){3732 if(control_type[i]==MaterialsRheologyBbarEnum || control_type[i]==MaterialsRheologyDbarEnum){ 3734 3733 input=(Input*)material->inputs->GetInput(control_type[i]); _assert_(input); 3735 3734 } … … 3758 3757 Input* input=NULL; 3759 3758 3760 if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheology ZbarEnum){3759 if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyDbarEnum){ 3761 3760 input=(Input*)material->inputs->GetInput(enum_type); 3762 3761 } … … 3776 3775 Input* input=NULL; 3777 3776 3778 if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheology ZbarEnum){3777 if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyDbarEnum){ 3779 3778 input=(Input*)material->inputs->GetInput(enum_type); 3780 3779 } … … 3795 3794 Input* input=NULL; 3796 3795 3797 if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheology ZbarEnum){3796 if(enum_type==MaterialsRheologyBbarEnum || enum_type==MaterialsRheologyDbarEnum){ 3798 3797 input=(Input*)material->inputs->GetInput(enum_type); 3799 3798 } … … 3826 3825 GradjBSSA(gradient,control_index); 3827 3826 break; 3828 case MaterialsRheology ZbarEnum:3829 Gradj ZSSA(gradient,control_index);3827 case MaterialsRheologyDbarEnum: 3828 GradjDSSA(gradient,control_index); 3830 3829 break; 3831 3830 case BalancethicknessThickeningRateEnum: … … 3919 3918 } 3920 3919 /*}}}*/ 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 /*}}}*/3961 3920 /*FUNCTION Tria::GradjBSSA{{{*/ 3962 3921 void Tria::GradjBSSA(Vector<IssmDouble>* gradient,int control_index){ … … 4016 3975 } 4017 3976 /*}}}*/ 4018 /*FUNCTION Tria::Gradj ZSSA{{{*/4019 void Tria::Gradj ZSSA(Vector<IssmDouble>* gradient,int control_index){3977 /*FUNCTION Tria::GradjDSSA{{{*/ 3978 void Tria::GradjDSSA(Vector<IssmDouble>* gradient,int control_index){ 4020 3979 4021 3980 /*Intermediaries*/ … … 4024 3983 IssmDouble vx,vy,lambda,mu,thickness,Jdet; 4025 3984 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]; 4027 3986 IssmDouble xyz_list[NUMVERTICES][3]; 4028 3987 IssmDouble basis[3],epsilon[3]; … … 4040 3999 Input* adjointx_input=inputs->GetInput(AdjointxEnum); _assert_(adjointx_input); 4041 4000 Input* adjointy_input=inputs->GetInput(AdjointyEnum); _assert_(adjointy_input); 4042 Input* rheology z_input=material->inputs->GetInput(MaterialsRheologyZbarEnum); _assert_(rheologyz_input);4001 Input* rheologyd_input=material->inputs->GetInput(MaterialsRheologyDbarEnum); _assert_(rheologyd_input); 4043 4002 4044 4003 /* Start looping on the number of gaussian points: */ … … 4049 4008 4050 4009 thickness_input->GetInputValue(&thickness,gauss); 4051 rheologyz_input->GetInputDerivativeValue(&dZ[0],&xyz_list[0][0],gauss);4052 4010 vx_input->GetInputDerivativeValue(&dvx[0],&xyz_list[0][0],gauss); 4053 4011 vy_input->GetInputDerivativeValue(&dvy[0],&xyz_list[0][0],gauss); … … 4056 4014 4057 4015 this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input); 4058 material->GetViscosity ZComplement(&viscosity_complement,&epsilon[0]);4016 material->GetViscosityDComplement(&viscosity_complement,&epsilon[0]); 4059 4017 4060 4018 GetJacobianDeterminant(&Jdet, &xyz_list[0][0],gauss); … … 5421 5379 5422 5380 /*Get input (either in element or material)*/ 5423 if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheology ZbarEnum){5381 if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheologyDbarEnum){ 5424 5382 input=(Input*)material->inputs->GetInput(control_enum); _assert_(input); 5425 5383 } … … 5456 5414 new_input = new TriaInput(control_enum,values,P1Enum); 5457 5415 5458 if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheology ZbarEnum){5416 if(control_enum==MaterialsRheologyBbarEnum || control_enum==MaterialsRheologyDbarEnum){ 5459 5417 input=(Input*)material->inputs->GetInput(control_enum); _assert_(input); 5460 5418 } -
issm/trunk-jpl/src/c/classes/Elements/Tria.h
r16125 r16146 139 139 void Gradj(Vector<IssmDouble>* gradient,int control_type,int control_index); 140 140 void GradjBGradient(Vector<IssmDouble>* gradient,int weight_index,int control_index); 141 void Gradj ZGradient(Vector<IssmDouble>* gradient,int weight_index,int control_index);141 void GradjDGradient(Vector<IssmDouble>* gradient,int weight_index,int control_index); 142 142 void GradjBSSA(Vector<IssmDouble>* gradient,int control_index); 143 void Gradj ZSSA(Vector<IssmDouble>* gradient,int control_index);143 void GradjDSSA(Vector<IssmDouble>* gradient,int control_index); 144 144 void GradjDragSSA(Vector<IssmDouble>* gradient,int control_index); 145 145 void GradjDragFS(Vector<IssmDouble>* gradient,int control_index); -
issm/trunk-jpl/src/c/classes/Materials/Matdamageice.cpp
r16125 r16146 161 161 } 162 162 /*}}}*/ 163 /*FUNCTION Matdamageice::GetD {{{*/ 164 IssmDouble Matdamageice::GetD(){ 165 166 /*Output*/ 167 IssmDouble D; 168 169 inputs->GetInputAverage(&D,MaterialsRheologyDEnum); 170 return D; 171 } 172 /*}}}*/ 163 173 /*FUNCTION Matdamageice::GetZ {{{*/ 164 174 IssmDouble Matdamageice::GetZ(){ 165 175 166 176 /*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); 171 181 } 172 182 /*}}}*/ … … 175 185 176 186 /*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 {{{*/ 193 IssmDouble Matdamageice::GetDbar(){ 194 195 /*Output*/ 196 IssmDouble Dbar; 197 inputs->GetInputAverage(&Dbar,MaterialsRheologyDbarEnum); 198 return Dbar; 181 199 } 182 200 /*}}}*/ … … 216 234 void Matdamageice::GetViscosity2d(IssmDouble* pviscosity, IssmDouble* epsilon){ 217 235 /*From a string tensor and a material object, return viscosity, using Glen's flow law. 218 Z* B236 (1-D) * B 219 237 viscosity= ------------------------------------------------------------------- 220 238 2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n] … … 235 253 /*Intermediary: */ 236 254 IssmDouble A,e; 237 IssmDouble Btmp,B,n, Z;255 IssmDouble Btmp,B,n,D; 238 256 239 257 /*Get B and n*/ 240 258 Btmp=GetBbar(); 241 Z=GetZbar();259 D=GetDbar(); 242 260 n=GetN(); 243 B= Z*Btmp;261 B=(1-D)*Btmp; 244 262 245 263 if (n==1){ … … 284 302 /*Return viscosity accounting for steady state power law creep [Thomas and SSA, 1982]: 285 303 * 286 * 304 * (1-D)* B 287 305 * viscosity3d= ------------------------------------------------------------------- 288 306 * 2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n] … … 303 321 /*Intermediaries: */ 304 322 IssmDouble A,e; 305 IssmDouble B,n, Z;323 IssmDouble B,n,D; 306 324 307 325 /*Get B, Z and n*/ 308 326 n=GetN(); 309 Z=GetZ();310 B= Z*GetB();327 D=GetD(); 328 B=(1-D)*GetB(); 311 329 312 330 if (n==1){ … … 355 373 /*Return viscosity accounting for steady state power law creep [Thomas and SSA, 1982]: 356 374 * 357 * B375 * (1-D)* B 358 376 * viscosity3d= ------------------------------------------------------------------- 359 377 * 2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n] … … 374 392 /*Intermediaries: */ 375 393 IssmDouble A,e; 376 IssmDouble B,n, Z;394 IssmDouble B,n,D; 377 395 IssmDouble eps0; 378 396 … … 380 398 eps0=pow(10.,(IssmDouble)-27); 381 399 n=GetN(); 382 Z=GetZ();383 B= Z*GetB();400 D=GetD(); 401 B=(1-D)*GetB(); 384 402 385 403 if (n==1){ … … 480 498 } 481 499 /*}}}*/ 482 /*FUNCTION Matdamageice::GetViscosity ZComplement {{{*/483 void Matdamageice::GetViscosity ZComplement(IssmDouble* pviscosity_complement, IssmDouble* epsilon){484 /*Return viscosity derivative for control method d(mu)/d Z:500 /*FUNCTION Matdamageice::GetViscosityDComplement {{{*/ 501 void Matdamageice::GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* epsilon){ 502 /*Return viscosity derivative for control method d(mu)/dD: 485 503 * 486 504 * B 487 * dviscosity= - ------------------------------------------------------------------505 * dviscosity= - ------------------------------------------------------------------- 488 506 * 2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n] 489 507 * … … 515 533 if(A==0){ 516 534 /*Maximum viscosity_complement for 0 shear areas: */ 517 viscosity_complement= 2.25*pow(10.,17);535 viscosity_complement=- 2.25*pow(10.,17); 518 536 } 519 537 else{ 520 538 e=(n-1)/(2*n); 521 539 522 viscosity_complement= B/(2*pow(A,e));540 viscosity_complement=- B/(2*pow(A,e)); 523 541 } 524 542 } 525 543 else{ 526 viscosity_complement= 4.5*pow(10.,17);544 viscosity_complement=- 4.5*pow(10.,17); 527 545 } 528 546 … … 530 548 _assert_(B>0); 531 549 _assert_(n>0); 532 _assert_(viscosity_complement >0);550 _assert_(viscosity_complement<0); 533 551 534 552 /*Return: */ … … 779 797 } 780 798 781 /*Get Z*/782 if (iomodel->Data(MaterialsRheology ZEnum)) {783 for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheology ZEnum)[iomodel->elements[num_vertices*index+i]-1];784 this->inputs->AddInput(new TriaInput(MaterialsRheology ZbarEnum,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)); 785 803 } 786 804 … … 799 817 } 800 818 break; 801 case MaterialsRheology ZbarEnum:802 if (iomodel->Data(MaterialsRheology ZEnum)){803 _assert_(iomodel->Data(MaterialsRheology ZEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));804 for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(MaterialsRheology ZEnum)[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]; 805 823 for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->Data(InversionMinParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i]; 806 824 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(MaterialsRheology ZbarEnum,TriaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));825 this->inputs->AddInput(new ControlInput(MaterialsRheologyDbarEnum,TriaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1)); 808 826 } 809 827 break; … … 837 855 } 838 856 839 /*Get Z*/840 if (iomodel->Data(MaterialsRheology ZEnum)) {841 for(i=0;i<num_vertices;i++) nodeinputs[i]=iomodel->Data(MaterialsRheology ZEnum)[iomodel->elements[num_vertices*index+i]-1];842 this->inputs->AddInput(new PentaInput(MaterialsRheology ZEnum,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)); 843 861 } 844 862 … … 857 875 } 858 876 break; 859 case MaterialsRheology ZbarEnum:860 if (iomodel->Data(MaterialsRheology ZEnum)){861 _assert_(iomodel->Data(MaterialsRheology ZEnum));_assert_(iomodel->Data(InversionMinParametersEnum)); _assert_(iomodel->Data(InversionMaxParametersEnum));862 for(j=0;j<num_vertices;j++)nodeinputs[j]=iomodel->Data(MaterialsRheology ZEnum)[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]; 863 881 for(j=0;j<num_vertices;j++)cmmininputs[j]=iomodel->Data(InversionMinParametersEnum)[(iomodel->elements[num_vertices*index+j]-1)*num_control_type+i]; 864 882 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(MaterialsRheology ZEnum,PentaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1));883 this->inputs->AddInput(new ControlInput(MaterialsRheologyDEnum,PentaInputEnum,nodeinputs,cmmininputs,cmmaxinputs,i+1)); 866 884 } 867 885 break; … … 885 903 name==MaterialsRheologyBbarEnum || 886 904 name==MaterialsRheologyNEnum || 887 name==MaterialsRheology ZEnum ||888 name==MaterialsRheology ZbarEnum905 name==MaterialsRheologyDEnum || 906 name==MaterialsRheologyDbarEnum 889 907 ){ 890 908 return true; -
issm/trunk-jpl/src/c/classes/Materials/Matdamageice.h
r16125 r16146 50 50 void GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon); 51 51 void GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon); 52 void GetViscosity ZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon);52 void GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon); 53 53 void GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon); 54 54 void GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon); … … 59 59 IssmDouble GetN(); 60 60 IssmDouble GetZ(); 61 IssmDouble GetD(); 61 62 IssmDouble GetZbar(); 63 IssmDouble GetDbar(); 62 64 bool IsInput(int name); 63 65 /*}}}*/ -
issm/trunk-jpl/src/c/classes/Materials/Material.h
r15654 r16146 30 30 virtual void GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon)=0; 31 31 virtual void GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon)=0; 32 virtual void GetViscosity ZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon)=0;32 virtual void GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon)=0; 33 33 virtual void GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon)=0; 34 34 virtual void GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon)=0; … … 39 39 virtual IssmDouble GetN()=0; 40 40 virtual IssmDouble GetZ()=0; 41 virtual IssmDouble GetD()=0; 41 42 virtual IssmDouble GetZbar()=0; 43 virtual IssmDouble GetDbar()=0; 42 44 43 45 }; -
issm/trunk-jpl/src/c/classes/Materials/Matice.h
r16125 r16146 57 57 void GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon); 58 58 void GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon); 59 void GetViscosityZComplement(IssmDouble*, IssmDouble*){_error_("not supported");};59 void GetViscosityDComplement(IssmDouble*, IssmDouble*){_error_("not supported");}; 60 60 void GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon); 61 61 void GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon); … … 65 65 IssmDouble GetBbar(); 66 66 IssmDouble GetZ(){_error_("not supported");}; 67 IssmDouble GetD(){_error_("not supported");}; 67 68 IssmDouble GetZbar(){_error_("not supported");}; 69 IssmDouble GetDbar(){_error_("not supported");}; 68 70 IssmDouble GetN(); 69 71 bool IsInput(int name); -
issm/trunk-jpl/src/c/classes/Materials/Matpar.h
r16125 r16146 87 87 void GetViscosity3dFS(IssmDouble* pviscosity3d, IssmDouble* epsilon){_error_("not supported");}; 88 88 void GetViscosityComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon){_error_("not supported");}; 89 void GetViscosity ZComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon){_error_("not supported");};89 void GetViscosityDComplement(IssmDouble* pviscosity_complement, IssmDouble* pepsilon){_error_("not supported");}; 90 90 void GetViscosityDerivativeEpsSquare(IssmDouble* pmu_prime, IssmDouble* pepsilon){_error_("not supported");}; 91 91 void GetViscosityDerivativeEpsSquareFS(IssmDouble* pmu_prime, IssmDouble* pepsilon){_error_("not supported");}; … … 95 95 IssmDouble GetBbar(){_error_("not supported");}; 96 96 IssmDouble GetN(){_error_("not supported");}; 97 IssmDouble GetD(){_error_("not supported");}; 97 98 IssmDouble GetZ(){_error_("not supported");}; 98 99 IssmDouble GetZbar(){_error_("not supported");}; 100 IssmDouble GetDbar(){_error_("not supported");}; 99 101 /*}}}*/ 100 102 /*Numerics: {{{*/ -
issm/trunk-jpl/src/c/modules/ModelProcessorx/Control/UpdateElementsAndMaterialsControl.cpp
r16050 r16146 40 40 case FrictionCoefficientEnum: iomodel->FetchData(1,FrictionCoefficientEnum); break; 41 41 case MaterialsRheologyBbarEnum: iomodel->FetchData(1,MaterialsRheologyBEnum); break; 42 case MaterialsRheology ZbarEnum: iomodel->FetchData(1,MaterialsRheologyZEnum); break;42 case MaterialsRheologyDbarEnum: iomodel->FetchData(1,MaterialsRheologyDEnum); break; 43 43 default: _error_("Control " << EnumToStringx(reCast<int,IssmDouble>(iomodel->Data(InversionControlParametersEnum)[i])) << " not implemented yet"); 44 44 } … … 59 59 60 60 /*Free data: */ 61 iomodel->DeleteData(4+7,InversionControlParametersEnum,InversionCostFunctionsCoefficientsEnum,InversionMinParametersEnum,InversionMaxParametersEnum,BalancethicknessThickeningRateEnum,VxEnum,VyEnum,ThicknessEnum,FrictionCoefficientEnum,MaterialsRheologyBEnum,MaterialsRheology ZEnum);61 iomodel->DeleteData(4+7,InversionControlParametersEnum,InversionCostFunctionsCoefficientsEnum,InversionMinParametersEnum,InversionMaxParametersEnum,BalancethicknessThickeningRateEnum,VxEnum,VyEnum,ThicknessEnum,FrictionCoefficientEnum,MaterialsRheologyBEnum,MaterialsRheologyDEnum); 62 62 } -
issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateElementsVerticesAndMaterials.cpp
r15941 r16146 64 64 break; 65 65 case MatdamageiceEnum: 66 iomodel->FetchData(3,MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheology ZEnum);66 iomodel->FetchData(3,MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheologyDEnum); 67 67 for (i=0;i<iomodel->numberofelements;i++) if(iomodel->my_elements[i]) materials->AddObject(new Matdamageice(i+1,i,iomodel)); 68 68 break; … … 73 73 /*Free data: */ 74 74 iomodel->DeleteData(8,MeshUpperelementsEnum,MeshLowerelementsEnum, 75 MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheology ZEnum,InversionControlParametersEnum,InversionMinParametersEnum,75 MaterialsRheologyBEnum,MaterialsRheologyNEnum,MaterialsRheologyDEnum,InversionControlParametersEnum,InversionMinParametersEnum, 76 76 InversionMaxParametersEnum); 77 77 -
issm/trunk-jpl/src/c/modules/ModelProcessorx/Stressbalance/UpdateElementsStressbalance.cpp
r16144 r16146 102 102 iomodel->FetchDataToInput(elements,LoadingforceYEnum); 103 103 if(materials_type==MatdamageiceEnum){ 104 iomodel->FetchDataToInput(elements,MaterialsRheology ZEnum);104 iomodel->FetchDataToInput(elements,MaterialsRheologyDEnum); 105 105 } 106 106 if(iomodel->dim==3){ -
issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
r16142 r16146 154 154 MaterialsRheologyZEnum, 155 155 MaterialsRheologyZbarEnum, 156 MaterialsRheologyDEnum, 157 MaterialsRheologyDbarEnum, 156 158 MaterialsRhoIceEnum, 157 159 MaterialsRhoWaterEnum, -
issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp
r16142 r16146 162 162 case MaterialsRheologyZEnum : return "MaterialsRheologyZ"; 163 163 case MaterialsRheologyZbarEnum : return "MaterialsRheologyZbar"; 164 case MaterialsRheologyDEnum : return "MaterialsRheologyD"; 165 case MaterialsRheologyDbarEnum : return "MaterialsRheologyDbar"; 164 166 case MaterialsRhoIceEnum : return "MaterialsRhoIce"; 165 167 case MaterialsRhoWaterEnum : return "MaterialsRhoWater"; -
issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp
r16142 r16146 165 165 else if (strcmp(name,"MaterialsRheologyZ")==0) return MaterialsRheologyZEnum; 166 166 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; 167 169 else if (strcmp(name,"MaterialsRhoIce")==0) return MaterialsRhoIceEnum; 168 170 else if (strcmp(name,"MaterialsRhoWater")==0) return MaterialsRhoWaterEnum; … … 258 260 else if (strcmp(name,"GiaLithosphereThickness")==0) return GiaLithosphereThicknessEnum; 259 261 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;262 262 else stage=3; 263 263 } 264 264 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; 266 268 else if (strcmp(name,"TimesteppingTimeAdapt")==0) return TimesteppingTimeAdaptEnum; 267 269 else if (strcmp(name,"TimesteppingTimeStep")==0) return TimesteppingTimeStepEnum; … … 381 383 else if (strcmp(name,"Penpair")==0) return PenpairEnum; 382 384 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;385 385 else stage=4; 386 386 } 387 387 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; 389 391 else if (strcmp(name,"NodeSId")==0) return NodeSIdEnum; 390 392 else if (strcmp(name,"VectorParam")==0) return VectorParamEnum; … … 504 506 else if (strcmp(name,"P0")==0) return P0Enum; 505 507 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;508 508 else stage=5; 509 509 } 510 510 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; 512 514 else if (strcmp(name,"P2")==0) return P2Enum; 513 515 else if (strcmp(name,"P2xP1")==0) return P2xP1Enum; -
issm/trunk-jpl/src/m/classes/inversion.m
r15860 r16146 116 116 md = checkfield(md,'inversion.incomplete_adjoint','values',[0 1]); 117 117 md = checkfield(md,'inversion.control_parameters','cell',1,'values',... 118 {'BalancethicknessThickeningRate' 'FrictionCoefficient' 'MaterialsRheologyBbar' 'MaterialsRheology Zbar' 'Vx' 'Vy' 'Thickness'});118 {'BalancethicknessThickeningRate' 'FrictionCoefficient' 'MaterialsRheologyBbar' 'MaterialsRheologyDbar' 'Vx' 'Vy' 'Thickness'}); 119 119 md = checkfield(md,'inversion.nsteps','numel',1,'>=',1); 120 120 md = checkfield(md,'inversion.maxiter_per_step','size',[md.inversion.nsteps 1],'>=',0); -
issm/trunk-jpl/src/m/classes/inversion.py
r15860 r16146 124 124 md = checkfield(md,'inversion.tao','values',[0,1]) 125 125 md = checkfield(md,'inversion.incomplete_adjoint','values',[0,1]) 126 md = checkfield(md,'inversion.control_parameters','cell',1,'values',['BalancethicknessThickeningRate','FrictionCoefficient','MaterialsRheologyBbar','MaterialsRheology Zbar','Vx','Vy'])126 md = checkfield(md,'inversion.control_parameters','cell',1,'values',['BalancethicknessThickeningRate','FrictionCoefficient','MaterialsRheologyBbar','MaterialsRheologyDbar','Vx','Vy']) 127 127 md = checkfield(md,'inversion.nsteps','numel',[1],'>=',1) 128 128 md = checkfield(md,'inversion.maxiter_per_step','size',[md.inversion.nsteps],'>=',0) -
issm/trunk-jpl/src/m/classes/matdamageice.m
r15825 r16146 19 19 rheology_B = NaN; 20 20 rheology_n = NaN; 21 rheology_ Z= NaN;21 rheology_D = NaN; 22 22 rheology_law = ''; 23 23 … … 101 101 md = checkfield(md,'materials.rheology_B','>',0,'size',[md.mesh.numberofvertices 1]); 102 102 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]); 104 104 md = checkfield(md,'materials.rheology_law','values',{'None' 'Paterson' 'Arrhenius' 'LliboutryDuval'}); 105 105 … … 128 128 fielddisplay(obj,'rheology_B','flow law parameter [Pa/s^(1/n)]'); 129 129 fielddisplay(obj,'rheology_n','Glen''s flow law exponent'); 130 fielddisplay(obj,'rheology_ Z','rheology multiplier');130 fielddisplay(obj,'rheology_D','damage tensor (scalar)'); 131 131 fielddisplay(obj,'rheology_law',['law for the temperature dependance of the rheology: ''None'', ''Paterson'', ''Arrhenius'' or ''LliboutryDuval''']); 132 132 fielddisplay(obj,'lithosphere_shear_modulus','Lithosphere shear modulus [Pa]'); … … 151 151 WriteData(fid,'object',obj,'class','materials','fieldname','rheology_B','format','DoubleMat','mattype',1); 152 152 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); 154 154 WriteData(fid,'data',StringToEnum(obj.rheology_law),'enum',MaterialsRheologyLawEnum(),'format','Integer'); 155 155 -
issm/trunk-jpl/src/m/classes/matdamageice.py
r15852 r16146 27 27 self.rheology_B = float('NaN') 28 28 self.rheology_n = float('NaN') 29 self.rheology_ Z= float('NaN')29 self.rheology_D = float('NaN') 30 30 self.rheology_law = '' 31 31 … … 58 58 s+="%s\n" % fielddisplay(self,"rheology_B","flow law parameter [Pa/s^(1/n)]") 59 59 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)") 61 61 s+="%s\n" % fielddisplay(self,"rheology_law","law for the temperature dependance of the rheology: 'None', 'Paterson', 'Arrhenius' or 'LliboutryDuval'") 62 62 s+="%s\n" % fielddisplay(self,"lithosphere_shear_modulus","Lithosphere shear modulus [Pa]") … … 119 119 md = checkfield(md,'materials.rheology_B','>',0,'size',[md.mesh.numberofvertices]) 120 120 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]) 122 122 md = checkfield(md,'materials.rheology_law','values',['None','Paterson','Arrhenius','LliboutryDuval']) 123 123 md = checkfield(md,'materials.lithosphere_shear_modulus','>',0,'numel',[1]); … … 143 143 WriteData(fid,'object',self,'class','materials','fieldname','rheology_B','format','DoubleMat','mattype',1) 144 144 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) 146 146 WriteData(fid,'data',StringToEnum(self.rheology_law)[0],'enum',MaterialsRheologyLawEnum(),'format','Integer') 147 147 -
issm/trunk-jpl/src/m/enum/EnumDefinitions.py
r16142 r16146 154 154 def MaterialsRheologyZEnum(): return StringToEnum("MaterialsRheologyZ")[0] 155 155 def MaterialsRheologyZbarEnum(): return StringToEnum("MaterialsRheologyZbar")[0] 156 def MaterialsRheologyDEnum(): return StringToEnum("MaterialsRheologyD")[0] 157 def MaterialsRheologyDbarEnum(): return StringToEnum("MaterialsRheologyDbar")[0] 156 158 def MaterialsRhoIceEnum(): return StringToEnum("MaterialsRhoIce")[0] 157 159 def MaterialsRhoWaterEnum(): return StringToEnum("MaterialsRhoWater")[0] -
issm/trunk-jpl/test/NightlyRun/test270.m
r15771 r16146 5 5 md.materials.rheology_B=paterson(md.initialization.temperature); 6 6 md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1); 7 md.materials.rheology_ Z=0.5*ones(md.mesh.numberofvertices,1);7 md.materials.rheology_D=0.5*ones(md.mesh.numberofvertices,1); 8 8 md=setflowequation(md,'SSA','all'); 9 9 md.cluster=generic('name',oshostname(),'np',3); -
issm/trunk-jpl/test/NightlyRun/test270.py
r15771 r16146 17 17 md.materials.rheology_B=paterson(md.initialization.temperature) 18 18 md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1)) 19 md.materials.rheology_ Z=0.5*numpy.ones((md.mesh.numberofvertices,1))19 md.materials.rheology_D=0.5*numpy.ones((md.mesh.numberofvertices,1)) 20 20 md=setflowequation(md,'SSA','all') 21 21 md.cluster=generic('name',oshostname(),'np',3) -
issm/trunk-jpl/test/NightlyRun/test272.m
r15771 r16146 5 5 md.materials.rheology_B=paterson(md.initialization.temperature); 6 6 md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1); 7 md.materials.rheology_ Z=0.5*ones(md.mesh.numberofvertices,1);7 md.materials.rheology_D=0.5*ones(md.mesh.numberofvertices,1); 8 8 md=setflowequation(md,'SSA','all'); 9 md.verbose=verbose('control',true); 9 10 10 11 %control parameters 11 12 md.inversion.iscontrol=1; 12 md.inversion.control_parameters={'MaterialsRheology Zbar'};13 md.inversion.min_parameters= 10^-13*ones(md.mesh.numberofvertices,1);14 md.inversion.max_parameters= ones(md.mesh.numberofvertices,1);13 md.inversion.control_parameters={'MaterialsRheologyDbar'}; 14 md.inversion.min_parameters=zeros(md.mesh.numberofvertices,1); 15 md.inversion.max_parameters=(1-10^-13)*ones(md.mesh.numberofvertices,1); 15 16 md.inversion.nsteps=2; 16 17 md.inversion.cost_functions=101.*ones(md.inversion.nsteps,1); … … 26 27 27 28 %Fields and tolerances to track changes 28 field_names ={'Gradient','Misfits','MaterialsRheology Zbar','Pressure','Vel','Vx','Vy'};29 field_names ={'Gradient','Misfits','MaterialsRheologyDbar','Pressure','Vel','Vx','Vy'}; 29 30 field_tolerances={1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12}; 30 31 field_values={... 31 32 (md.results.StressbalanceSolution.Gradient1),... 32 33 (md.results.StressbalanceSolution.J),... 33 (md.results.StressbalanceSolution.MaterialsRheology Zbar),...34 (md.results.StressbalanceSolution.MaterialsRheologyDbar),... 34 35 (md.results.StressbalanceSolution.Pressure),... 35 36 (md.results.StressbalanceSolution.Vel),... -
issm/trunk-jpl/test/NightlyRun/test272.py
r15771 r16146 17 17 md.materials.rheology_B=paterson(md.initialization.temperature) 18 18 md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1)) 19 md.materials.rheology_ Z=0.5*numpy.ones((md.mesh.numberofvertices,1))19 md.materials.rheology_D=0.5*numpy.ones((md.mesh.numberofvertices,1)) 20 20 md=setflowequation(md,'SSA','all') 21 21 22 22 #control parameters 23 23 md.inversion.iscontrol=1 24 md.inversion.control_parameters=['MaterialsRheology Zbar']24 md.inversion.control_parameters=['MaterialsRheologyDbar'] 25 25 md.inversion.min_parameters=10**-13*numpy.ones((md.mesh.numberofvertices,1)) 26 26 md.inversion.max_parameters=numpy.ones((md.mesh.numberofvertices,1)) … … 38 38 39 39 #Fields and tolerances to track changes 40 field_names =['Gradient','Misfits','MaterialsRheology Zbar','Pressure','Vel','Vx','Vy']40 field_names =['Gradient','Misfits','MaterialsRheologyDbar','Pressure','Vel','Vx','Vy'] 41 41 field_tolerances=[1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12,1e-12] 42 42 field_values=[\ 43 43 md.results.StressbalanceSolution.Gradient1,\ 44 44 md.results.StressbalanceSolution.J,\ 45 md.results.StressbalanceSolution.MaterialsRheology Zbar,\45 md.results.StressbalanceSolution.MaterialsRheologyDbar,\ 46 46 md.results.StressbalanceSolution.Pressure,\ 47 47 md.results.StressbalanceSolution.Vel,\ -
issm/trunk-jpl/test/NightlyRun/test274.m
r16075 r16146 6 6 md.materials.rheology_B=paterson(md.initialization.temperature); 7 7 md.materials.rheology_n=3.*ones(md.mesh.numberofelements,1); 8 md.materials.rheology_ Z=0.5*ones(md.mesh.numberofvertices,1);8 md.materials.rheology_D=0.5*ones(md.mesh.numberofvertices,1); 9 9 md=setflowequation(md,'SSA','all'); 10 10 -
issm/trunk-jpl/test/NightlyRun/test274.py
r16075 r16146 19 19 md.materials.rheology_B=paterson(md.initialization.temperature) 20 20 md.materials.rheology_n=3.*numpy.ones((md.mesh.numberofelements,1)) 21 md.materials.rheology_ Z=0.5*numpy.ones((md.mesh.numberofvertices,1))21 md.materials.rheology_D=0.5*numpy.ones((md.mesh.numberofvertices,1)) 22 22 md=setflowequation(md,'SSA','all') 23 23
Note:
See TracChangeset
for help on using the changeset viewer.