Changeset 19091


Ignore:
Timestamp:
02/05/15 15:24:10 (10 years ago)
Author:
bdef
Message:

CHG: changes in how we deal with the zigzag counter on the EPL

Location:
issm/trunk-jpl/src/c
Files:
4 edited

Legend:

Unmodified
Added
Removed

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

    r19082 r19091  
    646646}/*}}}*/
    647647
    648 void  HydrologyDCEfficientAnalysis::HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask, int* eplzigzag_counter, Element* element){
     648void  HydrologyDCEfficientAnalysis::HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask, Vector<IssmDouble>* recurence, int* eplzigzag_counter, Element* element){
    649649
    650650        bool        active_element;
    … …  
    672672        int         eplflip_lock;
    673673        int         numnodes      =basalelement->GetNumberOfNodes();
     674        IssmDouble  new_active;
     675        IssmDouble  recurent;
    674676        IssmDouble* epl_thickness =xNew<IssmDouble>(numnodes);
    675677        IssmDouble* old_active    =xNew<IssmDouble>(numnodes);
    … …  
    697699
    698700        for(i=0;i<numnodes;i++){
    699                 /*Activate EPL if residual is >0 */
     701
     702                /*Dealing with the initial value to define zigzaging, preceding iteration the first time we see the node and current evolving vector after*/
     703                recurence->GetValue(&recurent,basalelement->nodes[i]->Sid());
     704                if (recurent>0)vec_mask->GetValue(&old_active[i],basalelement->nodes[i]->Sid());
     705                new_active=old_active[i];
     706                recurence->SetValue(basalelement->nodes[i]->Sid(),1.,ADD_VAL);
     707
     708                /*Now starting to look at the activations*/
    700709                if(residual[i]>0.){
    701                         if(old_active[i]==0) eplzigzag_counter[basalelement->nodes[i]->Lid()] ++;
    702710                        vec_mask->SetValue(basalelement->nodes[i]->Sid(),1.,INS_VAL);
    703                 }
    704 
     711                        new_active=1.0;
     712                }
    705713                /*If mask was already one, keep one or colapse*/
    706714                else if(old_active[i]>0.){
    707                        
    708715                        vec_mask->SetValue(basalelement->nodes[i]->Sid(),1.,INS_VAL);
     716                        new_active=1.0;
    709717                        /*If epl thickness gets under colapse thickness, close the layer*/
    710718                        if(epl_thickness[i]<colapse_thick){
    711                                 eplzigzag_counter[basalelement->nodes[i]->Lid()] ++;
    712                                 /*Avoid flipfloping between open and closed states*/
    713                                 if(eplzigzag_counter[basalelement->nodes[i]->Lid()]<eplflip_lock |eplflip_lock==0){
    714                                         vec_mask->SetValue(basalelement->nodes[i]->Sid(),0.,INS_VAL);
    715                                         epl_thickness[i]=init_thick;
    716                                 }
     719                                vec_mask->SetValue(basalelement->nodes[i]->Sid(),0.,INS_VAL);
     720                                new_active=0.0;
    717721                        }
     722                }
     723
     724                if(basalelement->nodes[i]->Sid()==42){
     725                        printf("old active is %e, and new active is %e, with counter= %i\n",old_active[i],new_active,eplzigzag_counter[basalelement->nodes[i]->Lid()]);
    718726                }
    719727                /*Increase of the efficient system is needed if the epl head reach the maximum value (sediment max value for now)*/
    … …  
    723731                                /*Increase of the domain is on the downstream node in term of sediment head*/
    724732                                if(sedhead[j] == sedheadmin){
    725                                         if(old_active[j]==0) {
    726                                                 eplzigzag_counter[basalelement->nodes[j]->Lid()] ++;
    727                                                 if(eplzigzag_counter[basalelement->nodes[i]->Lid()]<eplflip_lock |eplflip_lock==0){
    728                                                         vec_mask->SetValue(basalelement->nodes[j]->Sid(),1.,INS_VAL);
    729                                                 }
    730                                                 else{
    731                                                         vec_mask->SetValue(basalelement->nodes[j]->Sid(),0.,INS_VAL);
    732                                                 }
    733                                         }
     733                                        vec_mask->SetValue(basalelement->nodes[j]->Sid(),1.,INS_VAL);
    734734                                }
    735735                        }
     736                }
     737                /*Now checking for nodes zigzaging between open and close*/
     738                if (old_active[i] != new_active){
     739                        eplzigzag_counter[basalelement->nodes[i]->Lid()] ++;
     740                        /*If node changed too much of state, fix it to it's last known state*/
     741                        if(eplzigzag_counter[basalelement->nodes[i]->Lid()]>eplflip_lock & eplflip_lock!=0){
     742                                vec_mask->SetValue(basalelement->nodes[i]->Sid(),old_active[i],INS_VAL);
     743                                new_active=old_active[i];
     744                        }
     745                }
     746                        /*If node is now closed bring its thickness back to initial*/
     747                if (new_active==0.){
     748                        epl_thickness[i]=init_thick;
    736749                }
    737750        }

  • issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.h

    r18719 r19091  
    4040                IssmDouble GetHydrologyKMatrixTransfer(Element* element, Gauss* gauss, Input* sed_head_input, Input* epl_head_input, Input* thick_input, Input* base_input);
    4141                IssmDouble GetHydrologyPVectorTransfer(Element* element, Gauss* gauss, Input* sed_head_input, Input* epl_head_input, Input* thick_input, Input* base_input);
    42                 void HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask, int* eplzigzag_counter, Element* element);
     42                void HydrologyEPLGetMask(Vector<IssmDouble>* vec_mask, Vector<IssmDouble>* recurence, int* eplzigzag_counter, Element* element);
    4343                void HydrologyEPLGetActive(Vector<IssmDouble>* active_vec, Element* element);
    4444                void GetHydrologyDCInefficientHmax(IssmDouble* ph_max,Element* element, Node* innode);

  • issm/trunk-jpl/src/c/classes/FemModel.cpp

    r19087 r19091  
    19131913
    19141914        Vector<IssmDouble>* mask                                                        = NULL;
     1915        Vector<IssmDouble>* recurence                           = NULL;
    19151916        IssmDouble*         serial_mask                         = NULL;
    19161917        Vector<IssmDouble>* active                                              = NULL;
    … …  
    19221923        /*Step 1: update mask, the mask might be extended by residual and/or using downstream sediment head*/
    19231924        mask=new Vector<IssmDouble>(this->nodes->NumberOfNodes(HydrologyDCEfficientAnalysisEnum));
     1925        recurence=new Vector<IssmDouble>(this->nodes->NumberOfNodes(HydrologyDCEfficientAnalysisEnum));
    19241926        this->parameters->FindParam(&eplzigzag_counter,NULL,EplZigZagCounterEnum);
    1925 
    19261927        for (int i=0;i<elements->Size();i++){
    19271928                Element* element=xDynamicCast<Element*>(elements->GetObjectByOffset(i));
    1928                 effanalysis->HydrologyEPLGetMask(mask,eplzigzag_counter,element);
     1929                effanalysis->HydrologyEPLGetMask(mask,recurence,eplzigzag_counter,element);
    19291930        }
    19301931        this->parameters->SetParam(eplzigzag_counter,this->nodes->Size(),EplZigZagCounterEnum);
    … …  
    19341935        xDelete<int>(eplzigzag_counter);
    19351936        delete mask;
     1937        delete recurence;
    19361938
    19371939        /*Update Mask*/

  • issm/trunk-jpl/src/c/solutionsequences/solutionsequence_hydro_nonlinear.cpp

    r19009 r19091  
    9292                femmodel->parameters->SetParam(HydrologySedimentEnum,HydrologyLayerEnum);
    9393               
    94                 /*Reset constraint on the ZigZag Lock, this thing doesn't work, it have to disapear*/
     94                /*Reset constraint on the ZigZag Lock*/
    9595                ResetConstraintsx(femmodel);
    9696
Note: See TracChangeset for help on using the changeset viewer.