Changeset 11229


Ignore:
Timestamp:
01/26/12 15:36:05 (13 years ago)
Author:
seroussi
Message:

added possibility to compute stress tensor components when specified in requested_outputs

Location:
issm/trunk-jpl/src
Files:
7 added
7 edited

Legend:

Unmodified
Added
Removed

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

    r11199 r11229  
    366366        HydrologyWaterVxEnum,
    367367        HydrologyWaterVyEnum,
     368        StressTensorEnum,
     369        StressTensorxxEnum,
     370        StressTensorxyEnum,
     371        StressTensorxzEnum,
     372        StressTensoryyEnum,
     373        StressTensoryzEnum,
     374        StressTensorzzEnum,
    368375        IceVolumeEnum,
    369376        /*}}}*/

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

    r11202 r11229  
    356356                case HydrologyWaterVxEnum : return "HydrologyWaterVx";
    357357                case HydrologyWaterVyEnum : return "HydrologyWaterVy";
     358                case StressTensorEnum : return "StressTensor";
     359                case StressTensorxxEnum : return "StressTensorxx";
     360                case StressTensorxyEnum : return "StressTensorxy";
     361                case StressTensorxzEnum : return "StressTensorxz";
     362                case StressTensoryyEnum : return "StressTensoryy";
     363                case StressTensoryzEnum : return "StressTensoryz";
     364                case StressTensorzzEnum : return "StressTensorzz";
    358365                case IceVolumeEnum : return "IceVolume";
    359366                case P0Enum : return "P0";

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

    r11202 r11229  
    354354        else if (strcmp(name,"HydrologyWaterVx")==0) return HydrologyWaterVxEnum;
    355355        else if (strcmp(name,"HydrologyWaterVy")==0) return HydrologyWaterVyEnum;
     356        else if (strcmp(name,"StressTensor")==0) return StressTensorEnum;
     357        else if (strcmp(name,"StressTensorxx")==0) return StressTensorxxEnum;
     358        else if (strcmp(name,"StressTensorxy")==0) return StressTensorxyEnum;
     359        else if (strcmp(name,"StressTensorxz")==0) return StressTensorxzEnum;
     360        else if (strcmp(name,"StressTensoryy")==0) return StressTensoryyEnum;
     361        else if (strcmp(name,"StressTensoryz")==0) return StressTensoryzEnum;
     362        else if (strcmp(name,"StressTensorzz")==0) return StressTensorzzEnum;
    356363        else if (strcmp(name,"IceVolume")==0) return IceVolumeEnum;
    357364        else if (strcmp(name,"P0")==0) return P0Enum;

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

    r11192 r11229  
    477477        _error_("Not implemented yet");
    478478
     479}
     480/*}}}*/
     481/*FUNCTION Penta::ComputeStressTensor {{{1*/
     482void  Penta::ComputeStressTensor(){
     483
     484        int         iv;
     485        double      xyz_list[NUMVERTICES][3];
     486        double      pressure,viscosity;
     487        double      epsilon[6]; /* epsilon=[exx,eyy,exy];*/
     488        double      sigma_xx[NUMVERTICES];
     489        double          sigma_yy[NUMVERTICES];
     490        double          sigma_zz[NUMVERTICES];
     491        double      sigma_xy[NUMVERTICES];
     492        double          sigma_xz[NUMVERTICES];
     493        double          sigma_yz[NUMVERTICES];
     494        GaussPenta* gauss=NULL;
     495
     496        /*retrive parameters: */
     497        parameters->FindParam(&analysis_type,AnalysisTypeEnum);
     498
     499        /* Get node coordinates and dof list: */
     500        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
     501
     502        /*Retrieve all inputs we will be needing: */
     503        Input* pressure_input=inputs->GetInput(PressureEnum); _assert_(pressure_input);
     504        Input* vx_input=inputs->GetInput(VxEnum);             _assert_(vx_input);
     505        Input* vy_input=inputs->GetInput(VyEnum);             _assert_(vy_input);
     506        Input* vz_input=inputs->GetInput(VzEnum);             _assert_(vz_input);
     507
     508        /* Start looping on the number of vertices: */
     509        gauss=new GaussPenta();
     510        for (int iv=0;iv<NUMVERTICES;iv++){
     511                gauss->GaussVertex(iv);
     512
     513                /*Compute strain rate viscosity and pressure: */
     514                this->GetStrainRate3d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input,vz_input);
     515                matice->GetViscosity3d(&viscosity,&epsilon[0]);
     516                pressure_input->GetInputValue(&pressure,gauss);
     517
     518                /*Compute Stress*/
     519                sigma_xx[iv]=2*viscosity*epsilon[0]-pressure; // sigma = nu eps - pressure
     520                sigma_yy[iv]=2*viscosity*epsilon[1]-pressure;
     521                sigma_zz[iv]=2*viscosity*epsilon[2]-pressure;
     522                sigma_xy[iv]=2*viscosity*epsilon[3];
     523                sigma_xz[iv]=2*viscosity*epsilon[4];
     524                sigma_yz[iv]=2*viscosity*epsilon[5];
     525        }
     526       
     527        /*Add Stress tensor components into inputs*/
     528        this->inputs->AddInput(new PentaVertexInput(StressTensorxxEnum,&sigma_xx[0]));
     529        this->inputs->AddInput(new PentaVertexInput(StressTensorxyEnum,&sigma_xy[0]));
     530        this->inputs->AddInput(new PentaVertexInput(StressTensorxzEnum,&sigma_xz[0]));
     531        this->inputs->AddInput(new PentaVertexInput(StressTensoryyEnum,&sigma_yy[0]));
     532        this->inputs->AddInput(new PentaVertexInput(StressTensoryzEnum,&sigma_yz[0]));
     533        this->inputs->AddInput(new PentaVertexInput(StressTensorzzEnum,&sigma_zz[0]));
     534
     535        /*Clean up and return*/
     536        delete gauss;
    479537}
    480538/*}}}*/
    … …  
    24252483                                /*erase input: */
    24262484                                inputs->DeleteInput(output_enum);
     2485                                break;
     2486
     2487                        case StressTensorEnum:
     2488                                this->ComputeStressTensor();
     2489                                InputToResult(StressTensorxxEnum,step,time);
     2490                                InputToResult(StressTensorxyEnum,step,time);
     2491                                InputToResult(StressTensorxzEnum,step,time);
     2492                                InputToResult(StressTensoryyEnum,step,time);
     2493                                InputToResult(StressTensoryzEnum,step,time);
     2494                                InputToResult(StressTensorzzEnum,step,time);
    24272495                                break;
    24282496

  • issm/trunk-jpl/src/c/objects/Elements/Penta.h

    r11001 r11229  
    8383                void   ComputeBasalStress(Vec sigma_b);
    8484                void   ComputeStrainRate(Vec eps);
     85                void   ComputeStressTensor();
    8586                void   Configure(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters);
    8687                void   SetCurrentConfiguration(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters);

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

    r11192 r11229  
    896896void  Tria::ComputeStrainRate(Vec eps){
    897897        _error_("Not Implemented yet");
     898}
     899/*}}}*/
     900/*FUNCTION Tria::ComputeStressTensor {{{1*/
     901void  Tria::ComputeStressTensor(){
     902
     903        int         iv;
     904        double      xyz_list[NUMVERTICES][3];
     905        double      pressure,viscosity;
     906        double      epsilon[3]; /* epsilon=[exx,eyy,exy];*/
     907        double      sigma_xx[NUMVERTICES];
     908        double          sigma_yy[NUMVERTICES];
     909        double          sigma_zz[NUMVERTICES]={0,0,0};
     910        double      sigma_xy[NUMVERTICES];
     911        double          sigma_xz[NUMVERTICES]={0,0,0};
     912        double          sigma_yz[NUMVERTICES]={0,0,0};
     913        GaussTria* gauss=NULL;
     914
     915        /* Get node coordinates and dof list: */
     916        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
     917
     918        /*Retrieve all inputs we will be needing: */
     919        Input* pressure_input=inputs->GetInput(PressureEnum); _assert_(pressure_input);
     920        Input* vx_input=inputs->GetInput(VxEnum);             _assert_(vx_input);
     921        Input* vy_input=inputs->GetInput(VyEnum);             _assert_(vy_input);
     922
     923        /* Start looping on the number of vertices: */
     924        gauss=new GaussTria();
     925        for (int iv=0;iv<NUMVERTICES;iv++){
     926                gauss->GaussVertex(iv);
     927
     928                /*Compute strain rate viscosity and pressure: */
     929                this->GetStrainRate2d(&epsilon[0],&xyz_list[0][0],gauss,vx_input,vy_input);
     930                matice->GetViscosity2d(&viscosity,&epsilon[0]);
     931                pressure_input->GetInputValue(&pressure,gauss);
     932
     933                /*Compute Stress*/
     934                sigma_xx[iv]=2*viscosity*epsilon[0]-pressure; // sigma = nu eps - pressure
     935                sigma_yy[iv]=2*viscosity*epsilon[1]-pressure;
     936                sigma_xy[iv]=2*viscosity*epsilon[2];
     937        }
     938       
     939        /*Add Stress tensor components into inputs*/
     940        this->inputs->AddInput(new TriaVertexInput(StressTensorxxEnum,&sigma_xx[0]));
     941        this->inputs->AddInput(new TriaVertexInput(StressTensorxyEnum,&sigma_xy[0]));
     942        this->inputs->AddInput(new TriaVertexInput(StressTensorxzEnum,&sigma_xz[0]));
     943        this->inputs->AddInput(new TriaVertexInput(StressTensoryyEnum,&sigma_yy[0]));
     944        this->inputs->AddInput(new TriaVertexInput(StressTensoryzEnum,&sigma_yz[0]));
     945        this->inputs->AddInput(new TriaVertexInput(StressTensorzzEnum,&sigma_zz[0]));
     946
     947        /*Clean up and return*/
     948        delete gauss;
    898949}
    899950/*}}}*/
    … …  
    21342185                /*this input does not exist, compute it, and then transfer to results: */
    21352186                switch(output_enum){
     2187                        case StressTensorEnum:
     2188                                this->ComputeStressTensor();
     2189                                InputToResult(StressTensorxxEnum,step,time);
     2190                                InputToResult(StressTensorxyEnum,step,time);
     2191                                InputToResult(StressTensorxzEnum,step,time);
     2192                                InputToResult(StressTensoryyEnum,step,time);
     2193                                InputToResult(StressTensoryzEnum,step,time);
     2194                                InputToResult(StressTensorzzEnum,step,time);
     2195                                break;
     2196
    21362197                        default:
    21372198                                /*do nothing, no need to derail the computation because one of the outputs requested cannot be found: */

  • issm/trunk-jpl/src/c/objects/Elements/Tria.h

    r11001 r11229  
    7979                void   ComputeBasalStress(Vec sigma_b);
    8080                void   ComputeStrainRate(Vec eps);
     81                void   ComputeStressTensor();
    8182                void   Configure(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters);
    8283                void   SetCurrentConfiguration(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters);
Note: See TracChangeset for help on using the changeset viewer.