Changeset 17434


Ignore:
Timestamp:
03/14/14 12:23:22 (11 years ago)
Author:
jbondzio
Message:

ADD: Coupling LSM to thermal and enthalpy analysis

Location:
issm/trunk-jpl/src/c/analyses
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp

    r17275 r17434  
    2525void EnthalpyAnalysis::UpdateElements(Elements* elements,IoModel* iomodel,int analysis_counter,int analysis_type){/*{{{*/
    2626
    27         bool dakota_analysis;
     27        bool dakota_analysis, islevelset;
    2828        bool isenthalpy;
    2929
    … …  
    4949
    5050        iomodel->Constant(&dakota_analysis,QmuIsdakotaEnum);
     51        iomodel->Constant(&islevelset,TransientIslevelsetEnum);
    5152
    5253        iomodel->FetchDataToInput(elements,ThicknessEnum);
    … …  
    8384                elements->InputDuplicate(VxMeshEnum,QmuVzMeshEnum);
    8485        }
     86        if(islevelset){
     87                iomodel->FetchDataToInput(elements,IceMaskNodeActivationEnum);
     88                iomodel->FetchDataToInput(elements,MeshVertexonbedEnum); // required for updating active nodes
     89        }
    8590
    8691        /*Free data: */
    … …  
    196201ElementMatrix* EnthalpyAnalysis::CreateKMatrix(Element* element){/*{{{*/
    197202
     203        /* Check if ice in element */
     204        if(!element->IsIceInElement()) return NULL;
     205
    198206        /*compute all stiffness matrices for this element*/
    199207        ElementMatrix* Ke1=CreateKMatrixVolume(element);
    … …  
    207215}/*}}}*/
    208216ElementMatrix* EnthalpyAnalysis::CreateKMatrixVolume(Element* element){/*{{{*/
     217
     218        /* Check if ice in element */
     219        if(!element->IsIceInElement()) return NULL;
    209220
    210221        /*Intermediaries */
    … …  
    339350ElementMatrix* EnthalpyAnalysis::CreateKMatrixShelf(Element* element){/*{{{*/
    340351
     352        /* Check if ice in element */
     353        if(!element->IsIceInElement()) return NULL;
     354
    341355        /*Initialize Element matrix and return if necessary*/
    342356        if(!element->IsOnBed() || !element->IsFloating()) return NULL;
    … …  
    388402ElementVector* EnthalpyAnalysis::CreatePVector(Element* element){/*{{{*/
    389403
     404        /* Check if ice in element */
     405        if(!element->IsIceInElement()) return NULL;
     406
    390407        /*compute all load vectors for this element*/
    391408        ElementVector* pe1=CreatePVectorVolume(element);
    … …  
    401418}/*}}}*/
    402419ElementVector* EnthalpyAnalysis::CreatePVectorVolume(Element* element){/*{{{*/
     420
     421        /* Check if ice in element */
     422        if(!element->IsIceInElement()) return NULL;
    403423
    404424        /*Intermediaries*/
    … …  
    515535ElementVector* EnthalpyAnalysis::CreatePVectorSheet(Element* element){/*{{{*/
    516536
     537        /* Check if ice in element */
     538        if(!element->IsIceInElement()) return NULL;
     539
    517540        /* implementation of the basal condition decision chart of Aschwanden 2012, Fig.5 */
    518541        if(!element->IsOnBed() || element->IsFloating()) return NULL;
    … …  
    601624}/*}}}*/
    602625ElementVector* EnthalpyAnalysis::CreatePVectorShelf(Element* element){/*{{{*/
     626
     627        /* Check if ice in element */
     628        if(!element->IsIceInElement()) return NULL;
    603629
    604630        /*Get basal element*/
    … …  
    819845}/*}}}*/
    820846void EnthalpyAnalysis::UpdateConstraints(FemModel* femmodel){/*{{{*/
    821         /*Default, do nothing*/
     847
     848        bool islevelset;
     849        femmodel->parameters->FindParam(&islevelset,TransientIslevelsetEnum);
     850        if(islevelset){
     851                _printf0_("   Updating active and non-active nodes for ThermalAnalysis \n");
     852                SetActiveNodesLSMx(femmodel);
     853        }
    822854        return;
    823855}/*}}}*/
    … …  
    829861        int solution_type, i;
    830862        bool computebasalmeltingrates=true;
    831         bool isdrainage=true;
     863        bool isdrainage=false;
    832864        bool updatebasalconstraints=true;
    833865
    … …  
    862894        /*Calculate the basal melt rates of the enthalpy model after Aschwanden 2012*/
    863895        /* melting rate is positive when melting, negative when refreezing*/
     896
     897        /* Check if ice in element */
     898        if(!element->IsIceInElement()) return;
     899
     900        /* Only compute melt rates at the base of grounded ice*/
     901        if(!element->IsOnBed() || element->IsFloating()) return;
    864902
    865903        /* Intermediaries */
    … …  
    877915        IssmDouble *xyz_list_base = NULL;
    878916        int        *pairindices   = NULL;
    879 
    880         /* Only compute melt rates at the base of grounded ice*/
    881         if(!element->IsOnBed() || element->IsFloating()) return;
    882917
    883918        /*Fetch parameters and inputs */
    … …  
    10061041void EnthalpyAnalysis::DrainWaterfractionIcecolumn(Element* element){/*{{{*/
    10071042
     1043        /* Check if ice in element */
     1044        if(!element->IsIceInElement()) return;
     1045
    10081046        /* Only drain waterfraction of ice column from element at base*/
    10091047        if(!element->IsOnBed()) return; //FIXME: allow freeze on for floating elements
    … …  
    10431081}/*}}}*/
    10441082void EnthalpyAnalysis::DrainWaterfraction(Element* element, IssmDouble* pdrainrate_element){/*{{{*/
     1083
     1084        /* Check if ice in element */
     1085        if(!element->IsIceInElement()) return;
    10451086
    10461087        /*Intermediaries*/
    … …  
    11021143void EnthalpyAnalysis::UpdateBasalConstraints(Element* element){/*{{{*/
    11031144
     1145        /* Check if ice in element */
     1146        if(!element->IsIceInElement()) return;
     1147
     1148        /* Only update Constraints at the base of grounded ice*/
     1149        if(!(element->IsOnBed()) || element->IsFloating()) return;
     1150
    11041151        /*Intermediary*/
    11051152        bool        isdynamicbasalspc,setspc;
    … …  
    11101157        int        *indices = NULL, *indicesup = NULL;
    11111158        Node*       node = NULL;
    1112        
    1113         /* Only update Constraints at the base of grounded ice*/
    1114         if(!(element->IsOnBed()) || element->IsFloating()) return;
    11151159
    11161160        /*Check wether dynamic basal boundary conditions are activated */

  • issm/trunk-jpl/src/c/analyses/ThermalAnalysis.cpp

    r17275 r17434  
    4343        }
    4444
    45         bool dakota_analysis;
     45        bool dakota_analysis, islevelset;
    4646        iomodel->Constant(&dakota_analysis,QmuIsdakotaEnum);
     47        iomodel->Constant(&islevelset,TransientIslevelsetEnum);
    4748
    4849        iomodel->FetchDataToInput(elements,ThicknessEnum);
    … …  
    7576                elements->InputDuplicate(VxMeshEnum,QmuVzMeshEnum);
    7677        }
     78        if(islevelset){
     79                iomodel->FetchDataToInput(elements,IceMaskNodeActivationEnum);
     80                iomodel->FetchDataToInput(elements,MeshVertexonbedEnum); // required for updating active nodes
     81        }
    7782}/*}}}*/
    7883void ThermalAnalysis::CreateNodes(Nodes* nodes,IoModel* iomodel){/*{{{*/
    … …  
    124129ElementMatrix* ThermalAnalysis::CreateKMatrix(Element* element){/*{{{*/
    125130
     131        /* Check if ice in element */
     132        if(!element->IsIceInElement()) return NULL;
     133
    126134        /*compute all stiffness matrices for this element*/
    127135        ElementMatrix* Ke1=CreateKMatrixVolume(element);
    … …  
    135143}/*}}}*/
    136144ElementMatrix* ThermalAnalysis::CreateKMatrixVolume(Element* element){/*{{{*/
     145
     146        /* Check if ice in element */
     147        if(!element->IsIceInElement()) return NULL;
    137148
    138149        /*Intermediaries */
    … …  
    266277ElementMatrix* ThermalAnalysis::CreateKMatrixShelf(Element* element){/*{{{*/
    267278
     279        /* Check if ice in element */
     280        if(!element->IsIceInElement()) return NULL;
     281
    268282        /*Initialize Element matrix and return if necessary*/
    269283        if(!element->IsOnBed() || !element->IsFloating()) return NULL;
    … …  
    316330}/*}}}*/
    317331ElementVector* ThermalAnalysis::CreatePVector(Element* element){/*{{{*/
     332       
     333        /* Check if ice in element */
     334        if(!element->IsIceInElement()) return NULL;
    318335
    319336        /*compute all load vectors for this element*/
    … …  
    330347}/*}}}*/
    331348ElementVector* ThermalAnalysis::CreatePVectorVolume(Element* element){/*{{{*/
     349
     350        /* Check if ice in element */
     351        if(!element->IsIceInElement()) return NULL;
    332352
    333353        /*Intermediaries*/
    … …  
    410430ElementVector* ThermalAnalysis::CreatePVectorSheet(Element* element){/*{{{*/
    411431
     432        /* Check if ice in element */
     433        if(!element->IsIceInElement()) return NULL;
     434
    412435        /* Geothermal flux on ice sheet base and basal friction */
    413436        if(!element->IsOnBed() || element->IsFloating()) return NULL;
    … …  
    468491}/*}}}*/
    469492ElementVector* ThermalAnalysis::CreatePVectorShelf(Element* element){/*{{{*/
     493
     494        /* Check if ice in element */
     495        if(!element->IsIceInElement()) return NULL;
    470496
    471497        IssmDouble  t_pmp,dt,Jdet,scalar_ocean,pressure;
    … …  
    683709}/*}}}*/
    684710void ThermalAnalysis::UpdateConstraints(FemModel* femmodel){/*{{{*/
    685         /*Default, do nothing*/
     711
     712        bool islevelset;
     713        femmodel->parameters->FindParam(&islevelset,TransientIslevelsetEnum);
     714        if(islevelset){
     715                _printf0_("   Updating active and non-active nodes for ThermalAnalysis \n");
     716                SetActiveNodesLSMx(femmodel);
     717        }
    686718        return;
    687719}/*}}}*/
Note: See TracChangeset for help on using the changeset viewer.