Changeset 16612


Ignore:
Timestamp:
11/04/13 16:07:12 (11 years ago)
Author:
nathanma
Message:

CHG: several changes regarding surface normal velocity computation at free surface for 2D vertical transient run

Location:
issm/trunk-jpl/src
Files:
10 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp

    r16542 r16612  
    196196                if(dakota_analysis)elements->InputDuplicate(VzEnum,QmuVzEnum);
    197197        }
     198        if(iomodel->meshtype==Mesh2DverticalEnum){
     199              iomodel->FetchDataToInput(elements,MeshVertexonsurfaceEnum);
     200        }
    198201        if(isFS){
    199202                iomodel->FetchDataToInput(elements,MeshVertexonbedEnum);
  • issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

    r16600 r16612  
    24552455                                input=this->inputs->GetInput(output_enum);
    24562456                                break;
     2457                case SurfaceNormalVelocityEnum:
     2458                              this->ComputeSurfaceNormalVelocity();
     2459                                input=this->inputs->GetInput(output_enum);
     2460                                break;
    24572461                        default:
    24582462                                _error_("input "<<EnumToStringx(output_enum)<<" not found in element");
     
    24862490
    24872491                        vector->SetValues(NUMVERTICES,&sidlist[0],&values[0],ADD_VAL);
     2492                       
     2493                        this->inputs->DeleteInput(SurfaceNormalVelocityEnum);
     2494
    24882495                        break;
    24892496                                        }
     
    26422649        *(surface_normal+1) = normal[1]/normal_norm;
    26432650        *(surface_normal+2) = normal[2]/normal_norm;
     2651}
     2652/*}}}*/
     2653/*FUNCTION Tria::ComputeSurfaceNormalVelocity{{{*/
     2654void Tria::ComputeSurfaceNormalVelocity(){
     2655
     2656  IssmDouble      sum,tangential_vector[2],normal_vector[2],time,ux,uy;
     2657  IssmDouble      normal_velocity[NUMVERTICES],xyz_list[NUMVERTICES][3];
     2658  IssmDouble      value[NUMVERTICES],verticesonsurface[NUMVERTICES];
     2659
     2660  for(int iv=0;iv<NUMVERTICES;iv++){
     2661    normal_velocity[iv]=0.;
     2662    value[iv]=0.;
     2663  }
     2664
     2665  GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES);
     2666
     2667  GaussTria* gauss=new GaussTria();
     2668  Input* slope_input=inputs->GetInput(SurfaceSlopeXEnum); _assert_(slope_input);
     2669  //  Input* slope_input= inputs->GetInput(SurfaceEnum); _assert_(slope_input);
     2670  Input* vx_input=inputs->GetInput(VxEnum); _assert_(vx_input);
     2671  Input* vy_input=inputs->GetInput(VyEnum); _assert_(vy_input);
     2672
     2673
     2674  /*Get list of nodes on surface*/
     2675  GetInputListOnVertices(&verticesonsurface[0],MeshVertexonsurfaceEnum); 
     2676  sum = verticesonsurface[0]+verticesonsurface[1]+verticesonsurface[2];
     2677  _assert_(sum==0. || sum==1. || sum==2.);
     2678
     2679  /*Compute normal velocity for surface nodes from L2 projected slope*/
     2680  for(int iv=0;iv<NUMVERTICES;iv++){
     2681    if(verticesonsurface[iv] == 1){
     2682        gauss->GaussNode(P1Enum,iv);   
     2683        slope_input->GetInputValue(&value[iv],gauss);
     2684        //      slope_input->GetInputDerivativeValue(&value[iv],&xyz_list[0][0],gauss);
     2685        vx_input->GetInputValue(&ux,gauss);
     2686        vy_input->GetInputValue(&uy,gauss);
     2687        tangential_vector[0]=sqrt(1./(pow(value[iv],2.)+1.));
     2688        tangential_vector[1]=value[iv]*tangential_vector[0];
     2689        normal_vector[0]=-1.*tangential_vector[1];
     2690        normal_vector[1]=tangential_vector[0];
     2691        normal_velocity[iv]=ux*normal_vector[0]+uy*normal_vector[1];
     2692    }
     2693  }
     2694
     2695  delete gauss;
     2696  this->inputs->AddInput(new TriaInput(SurfaceNormalVelocityEnum,&normal_velocity[0],P1Enum));
     2697
    26442698}
    26452699/*}}}*/
     
    32243278        Input* smb_input = inputs->GetInput(SurfaceforcingsMassBalanceEnum); _assert_(smb_input);
    32253279        smb_input->GetInputAverage(&smb);                                                                                                                                                                                               // average smb on element in m ice s-1
    3226    Total_Smb=rho_ice*base*smb;                                                                                                                                                                                                                  // smb on element in kg s-1
     3280        Total_Smb=rho_ice*base*smb;                                                                                                                                                                                                                     // smb on element in kg s-1
    32273281
    32283282        /*Return: */
  • issm/trunk-jpl/src/c/classes/Elements/Tria.h

    r16600 r16612  
    7070                void        ComputeStrainRate(Vector<IssmDouble>* eps);
    7171                void        ComputeStressTensor();
     72                void        ComputeSurfaceNormalVelocity();
    7273                void        Configure(Elements* elements,Loads* loads,Nodes* nodesin,Vertices* verticesin,Materials* materials,Parameters* parameters);
    7374                void        SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters);
     
    241242                Seg*             SpawnSeg(int index1,int index2);
    242243                void             SurfaceNormal(IssmDouble* surface_normal, IssmDouble xyz_list[3][3]);
    243 
    244244                #ifdef _HAVE_STRESSBALANCE_
    245245                ElementMatrix* CreateKMatrixStressbalanceSSA(void);
  • issm/trunk-jpl/src/c/cores/AnalysisConfiguration.cpp

    r16533 r16612  
    2626
    2727                case StressbalanceSolutionEnum:
    28                         numanalyses=4;
     28                        numanalyses=5;
    2929                        analyses=xNew<int>(numanalyses);
    3030                        analyses[0]=StressbalanceAnalysisEnum;
     
    3232                        analyses[2]=StressbalanceSIAAnalysisEnum;
    3333                        analyses[3]=L2ProjectionBaseAnalysisEnum;
     34                        analyses[4]=ExtrudeFromBaseAnalysisEnum;
    3435                        break;
    3536
  • issm/trunk-jpl/src/c/cores/stressbalance_core.cpp

    r16518 r16612  
    4848
    4949        /*Compute slopes: */
    50         if(isSIA) surfaceslope_core(femmodel);
     50        if(isSIA || (isFS && meshtype==Mesh2DverticalEnum)) surfaceslope_core(femmodel);
    5151        if(isFS){
    5252                bedslope_core(femmodel);
  • issm/trunk-jpl/src/c/cores/surfaceslope_core.cpp

    r16529 r16612  
    3232                solutionsequence_linear(femmodel);
    3333        }
     34        if(meshtype==Mesh2DverticalEnum){
     35              femmodel->parameters->SetParam(SurfaceSlopeXEnum,InputToExtrudeEnum);
     36                extrudefrombase_core(femmodel);
     37        }
    3438
    3539        if(save_results){
  • issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h

    r16600 r16612  
    321321        FreeSurfaceBaseAnalysisEnum,
    322322        FreeSurfaceTopAnalysisEnum,
     323        SurfaceNormalVelocityEnum,
    323324        ExtrudeFromBaseAnalysisEnum,
    324325        ExtrudeFromTopAnalysisEnum,
  • issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp

    r16600 r16612  
    327327                case FreeSurfaceBaseAnalysisEnum : return "FreeSurfaceBaseAnalysis";
    328328                case FreeSurfaceTopAnalysisEnum : return "FreeSurfaceTopAnalysis";
     329                case SurfaceNormalVelocityEnum : return "SurfaceNormalVelocity";
    329330                case ExtrudeFromBaseAnalysisEnum : return "ExtrudeFromBaseAnalysis";
    330331                case ExtrudeFromTopAnalysisEnum : return "ExtrudeFromTopAnalysis";
  • issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp

    r16600 r16612  
    333333              else if (strcmp(name,"FreeSurfaceBaseAnalysis")==0) return FreeSurfaceBaseAnalysisEnum;
    334334              else if (strcmp(name,"FreeSurfaceTopAnalysis")==0) return FreeSurfaceTopAnalysisEnum;
     335              else if (strcmp(name,"SurfaceNormalVelocity")==0) return SurfaceNormalVelocityEnum;
    335336              else if (strcmp(name,"ExtrudeFromBaseAnalysis")==0) return ExtrudeFromBaseAnalysisEnum;
    336337              else if (strcmp(name,"ExtrudeFromTopAnalysis")==0) return ExtrudeFromTopAnalysisEnum;
     
    382383              else if (strcmp(name,"Element")==0) return ElementEnum;
    383384              else if (strcmp(name,"ExternalResult")==0) return ExternalResultEnum;
    384               else if (strcmp(name,"FileParam")==0) return FileParamEnum;
    385385         else stage=4;
    386386   }
    387387   if(stage==4){
    388               if (strcmp(name,"Input")==0) return InputEnum;
     388              if (strcmp(name,"FileParam")==0) return FileParamEnum;
     389              else if (strcmp(name,"Input")==0) return InputEnum;
    389390              else if (strcmp(name,"IntInput")==0) return IntInputEnum;
    390391              else if (strcmp(name,"InputToExtrude")==0) return InputToExtrudeEnum;
     
    505506              else if (strcmp(name,"QmuTemperature")==0) return QmuTemperatureEnum;
    506507              else if (strcmp(name,"HydrologyWaterVx")==0) return HydrologyWaterVxEnum;
    507               else if (strcmp(name,"HydrologyWaterVy")==0) return HydrologyWaterVyEnum;
    508508         else stage=5;
    509509   }
    510510   if(stage==5){
    511               if (strcmp(name,"StressTensor")==0) return StressTensorEnum;
     511              if (strcmp(name,"HydrologyWaterVy")==0) return HydrologyWaterVyEnum;
     512              else if (strcmp(name,"StressTensor")==0) return StressTensorEnum;
    512513              else if (strcmp(name,"StressTensorxx")==0) return StressTensorxxEnum;
    513514              else if (strcmp(name,"StressTensorxy")==0) return StressTensorxyEnum;
  • issm/trunk-jpl/src/m/enum/EnumDefinitions.py

    r16600 r16612  
    319319def FreeSurfaceBaseAnalysisEnum(): return StringToEnum("FreeSurfaceBaseAnalysis")[0]
    320320def FreeSurfaceTopAnalysisEnum(): return StringToEnum("FreeSurfaceTopAnalysis")[0]
     321def SurfaceNormalVelocityEnum(): return StringToEnum("SurfaceNormalVelocity")[0]
    321322def ExtrudeFromBaseAnalysisEnum(): return StringToEnum("ExtrudeFromBaseAnalysis")[0]
    322323def ExtrudeFromTopAnalysisEnum(): return StringToEnum("ExtrudeFromTopAnalysis")[0]
Note: See TracChangeset for help on using the changeset viewer.