Changeset 16042


Ignore:
Timestamp:
08/30/13 11:44:43 (12 years ago)
Author:
Mathieu Morlighem
Message:

NEW: speeding up matrix preallocation process again, by keeping track of flags to revert

Location:
issm/trunk-jpl/src/c/classes
Files:
15 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/classes/Elements/Element.h

    r16026 r16042  
    3636                virtual void   Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters)=0;
    3737                virtual void   SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters)=0;
    38                 virtual void   SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int set1_enum,int set2_enum)=0;
     38                virtual void   SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum)=0;
    3939                virtual void   CreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>*  Kfs)=0;
    4040                virtual void   CreateDVector(Vector<IssmDouble>* df)=0;
    … …  
    4545                virtual int    GetNumberOfNodes(void)=0;
    4646                virtual void   GetNodesSidList(int* sidlist)=0;
     47                virtual void   GetNodesLidList(int* sidlist)=0;
    4748
    4849                virtual int    Sid()=0;

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

    r16037 r16042  
    13271327}
    13281328/*}}}*/
     1329/*FUNCTION Penta::GetNodesLidList{{{*/
     1330void Penta::GetNodesLidList(int* lidlist){
     1331
     1332        _assert_(lidlist);
     1333        _assert_(nodes);
     1334        int numnodes = this->NumberofNodes();
     1335
     1336        for(int i=0;i<numnodes;i++){
     1337                lidlist[i]=nodes[i]->Lid();
     1338        }
     1339}
     1340/*}}}*/
    13291341/*FUNCTION Penta::GetNumberOfNodes;{{{*/
    13301342int Penta::GetNumberOfNodes(void){
    … …  
    32513263/*}}}*/
    32523264/*FUNCTION Penta::SetwiseNodeConnectivity{{{*/
    3253 void Penta::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int set1_enum,int set2_enum){
     3265void Penta::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum){
    32543266
    32553267        /*Intermediaries*/
    … …  
    32673279                        /*flag current node so that no other element processes it*/
    32683280                        flags[this->nodes[i]->Lid()]=true;
     3281
     3282                        int counter=0;
     3283                        while(flagsindices[counter]>=0) counter++;
     3284                        flagsindices[counter]=this->nodes[i]->Lid();
    32693285
    32703286                        /*if node is clone, we have an off-diagonal non-zero, else it is a diagonal non-zero*/

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

    r16026 r16042  
    8080                void   Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    8181                void   SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters);
    82                 void   SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int set1_enum,int set2_enum);
     82                void   SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum);
    8383                void   CreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* Kfs);
    8484                void   CreateDVector(Vector<IssmDouble>* df);
    … …  
    8989                int    GetNodeIndex(Node* node);
    9090                void   GetNodesSidList(int* sidlist);
     91                void   GetNodesLidList(int* lidlist);
    9192                int    GetNumberOfNodes(void);
    9293                void   GetSolutionFromInputs(Vector<IssmDouble>* solution);

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

    r16036 r16042  
    119119/*Other*/
    120120/*FUNCTION Tria::SetwiseNodeConnectivity{{{*/
    121 void Tria::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int set1_enum,int set2_enum){
     121void Tria::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum){
    122122
    123123        /*Intermediaries*/
    … …  
    135135                        /*flag current node so that no other element processes it*/
    136136                        flags[this->nodes[i]->Lid()]=true;
     137
     138                        int counter=0;
     139                        while(flagsindices[counter]>=0) counter++;
     140                        flagsindices[counter]=this->nodes[i]->Lid();
    137141
    138142                        /*if node is clone, we have an off-diagonal non-zero, else it is a diagonal non-zero*/
    … …  
    11101114        for(int i=0;i<numnodes;i++){
    11111115                sidlist[i]=nodes[i]->Sid();
     1116        }
     1117}
     1118/*}}}*/
     1119/*FUNCTION Tria::GetNodesLidList{{{*/
     1120void Tria::GetNodesLidList(int* lidlist){
     1121
     1122        _assert_(lidlist);
     1123        _assert_(nodes);
     1124        int numnodes = this->NumberofNodes();
     1125
     1126        for(int i=0;i<numnodes;i++){
     1127                lidlist[i]=nodes[i]->Lid();
    11121128        }
    11131129}
    … …  
    75697585                        vx_input->GetInputValue(&vx,gauss);
    75707586                        vy_input->GetInputValue(&vy,gauss);
    7571                         norms = sqrt(slope[0]*slope[0]+slope[1]*slope[1]);
     7587                        norms = sqrt(slope[0]*slope[0]+slope[1]*slope[1]+1.e-10);
    75727588                        normv = sqrt(vx*vx + vy*vy);
    75737589                        if(normv>15./(365.*24.*3600.)) slope[0] = -vx/normv*norms;
    … …  
    76217637                        vx_input->GetInputValue(&vx,gauss);
    76227638                        vy_input->GetInputValue(&vy,gauss);
    7623                         norms = sqrt(slope[0]*slope[0]+slope[1]*slope[1]);
     7639                        norms = sqrt(slope[0]*slope[0]+slope[1]*slope[1]+1.e-10);
    76247640                        normv = sqrt(vx*vx + vy*vy);
    76257641                        if(normv>15./(365.*24.*3600.)) slope[1] = -vy/normv*norms;

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

    r16026 r16042  
    7676                void        Configure(Elements* elements,Loads* loads,Nodes* nodesin,Vertices* verticesin,Materials* materials,Parameters* parameters);
    7777                void        SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters);
    78                 void        SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int set1_enum,int set2_enum);
     78                void        SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum);
    7979                void        CreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* Kfs);
    8080                void        CreateDVector(Vector<IssmDouble>* df);
    … …  
    8484                int         GetNodeIndex(Node* node);
    8585                void        GetNodesSidList(int* sidlist);
     86                void        GetNodesLidList(int* lidlist);
    8687                int         GetNumberOfNodes(void);
    8788                int         Sid();

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

    r16038 r16042  
    518518        /*OK now count number of dofs and flag each nodes for each node i*/
    519519        bool *flags                  = xNew<bool>(localnumnodes);
     520        int  *flagsindices           = xNew<int>(localnumnodes);
    520521        int  *d_connectivity         = xNewZeroInit<int>(numnodes);
    521522        int  *o_connectivity         = xNewZeroInit<int>(numnodes);
    … …  
    523524        int  *all_connectivity_clone = xNewZeroInit<int>(numnodes);
    524525
     526        /*Resetting flags to false at eahc iteration takes a lot of time, so we keep track of the flags
     527         * to reset in flagsindices, initialized with -1*/
     528        for(i = 0;i<localnumnodes;i++) flags[i]        = false;
     529        for(i = 0;i<localnumnodes;i++) flagsindices[i] = -1;
     530
    525531        /*Create connectivity vector*/
    526532        for(i=0;i<nodes->Size();i++){
    … …  
    528534                if(node->InAnalysis(configuration_type)){
    529535
    530                         /*Reinitialize flags to 0*/
    531                         for(j=0;j<localnumnodes;j++) flags[j]=false;
     536                        /*Reinitialize flags to false*/
     537                        j=0;
     538                        while(true){
     539                                if(flagsindices[j]>=0){
     540                                        flags[flagsindices[j]] = false;
     541                                        flagsindices[j]        = -1;
     542                                        j++;
     543                                }
     544                                else{
     545                                        break;
     546                                }
     547                        }
     548
     549                        //for(j=0;j<localnumnodes;j++) flags[j]=false;
    532550
    533551                        /*Loop over elements that hold node number i*/
    534                         //if(head_e[node->Sid()]==-1 && head_l[node->Sid()]==-1){
    535                         //      printf("[%i] vertex %i\n",IssmComm::GetRank(),node->Sid()+1);
     552                        //if(head_e[node->Lid()]==-1 && head_l[node->Lid()]==-1){
     553                        //      printf("[%i] vertex %i\n",IssmComm::GetRank(),node->Lid()+1);
    536554                        //}
    537555                        for(j=head_e[node->Sid()];j!=-1;j=next_e[j]){
    538556                                offset=count2offset_e[j];
    539557                                element=dynamic_cast<Element*>(elements->GetObjectByOffset(offset));
    540                                 element->SetwiseNodeConnectivity(&d_nz,&o_nz,node,flags,set1enum,set2enum);
     558                                element->SetwiseNodeConnectivity(&d_nz,&o_nz,node,flags,flagsindices,set1enum,set2enum);
    541559                                if(node->IsClone()){
    542560                                        connectivity_clone[node->Sid()]+=d_nz+o_nz;
    … …  
    550568                                offset=count2offset_l[j];
    551569                                load=dynamic_cast<Load*>(loads->GetObjectByOffset(offset));
    552                                 load->SetwiseNodeConnectivity(&d_nz,&o_nz,node,flags,set1enum,set2enum);
     570                                load->SetwiseNodeConnectivity(&d_nz,&o_nz,node,flags,flagsindices,set1enum,set2enum);
    553571                                if(node->IsClone()){
    554572                                        connectivity_clone[node->Sid()]+=d_nz+o_nz;
    … …  
    562580        }
    563581        xDelete<bool>(flags);
     582        xDelete<int>(flagsindices);
    564583        xDelete<int>(count2offset_e);
    565584        xDelete<int>(head_e);

  • issm/trunk-jpl/src/c/classes/Loads/Load.h

    r15649 r16042  
    3030                virtual int   GetNumberOfNodes(void)=0;
    3131                virtual void  GetNodesSidList(int* sidlist)=0;
     32                virtual void  GetNodesLidList(int* lidlist)=0;
    3233                virtual void  SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters)=0;
    3334                virtual void  CreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* Kfs)=0;
    … …  
    3839                virtual void  PenaltyCreatePVector(Vector<IssmDouble>* pf, IssmDouble kmax)=0;
    3940                virtual bool  InAnalysis(int analysis_type)=0;
    40                 virtual void  SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int set1_enum,int set2_enum)=0;
     41                virtual void  SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum)=0;
    4142};
    4243#endif

  • issm/trunk-jpl/src/c/classes/Loads/Numericalflux.cpp

    r16035 r16042  
    325325}
    326326/*}}}*/
     327/*FUNCTION Numericalflux::GetNodesLidList{{{*/
     328void Numericalflux::GetNodesLidList(int* lidlist){
     329
     330        int type;
     331        inputs->GetInputValue(&type,NumericalfluxTypeEnum);
     332        _assert_(lidlist);
     333        _assert_(nodes);
     334
     335        switch(type){
     336                case InternalEnum:
     337                        for(int i=0;i<NUMNODES_INTERNAL;i++) lidlist[i]=nodes[i]->Lid();
     338                        return;
     339                case BoundaryEnum:
     340                        for(int i=0;i<NUMNODES_BOUNDARY;i++) lidlist[i]=nodes[i]->Lid();
     341                        return;
     342                default:
     343                        _error_("Numericalflux type " << EnumToStringx(type) << " not supported yet");
     344        }
     345}
     346/*}}}*/
    327347/*FUNCTION Numericalflux::GetNumberOfNodes{{{*/
    328348int Numericalflux::GetNumberOfNodes(void){
    … …  
    370390/*}}}*/
    371391/*FUNCTION Numericalflux::SetwiseNodeConnectivity{{{*/
    372 void Numericalflux::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int set1_enum,int set2_enum){
     392void Numericalflux::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum){
    373393
    374394        /*Output */
    … …  
    383403                        /*flag current node so that no other element processes it*/
    384404                        flags[this->nodes[i]->Lid()]=true;
     405
     406                        int counter=0;
     407                        while(flagsindices[counter]>=0) counter++;
     408                        flagsindices[counter]=this->nodes[i]->Lid();
    385409
    386410                        /*if node is clone, we have an off-diagonal non-zero, else it is a diagonal non-zero*/

  • issm/trunk-jpl/src/c/classes/Loads/Numericalflux.h

    r15767 r16042  
    6767                void CreatePVector(Vector<IssmDouble>* pf);
    6868                void GetNodesSidList(int* sidlist);
     69                void GetNodesLidList(int* lidlist);
    6970                int  GetNumberOfNodes(void);
    7071                void CreateJacobianMatrix(Matrix<IssmDouble>* Jff){_error_("Not implemented yet");};
    … …  
    7374                void PenaltyCreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* kfs, IssmDouble kmax);
    7475                void PenaltyCreatePVector(Vector<IssmDouble>* pf, IssmDouble kmax);
    75                 void SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int set1_enum,int set2_enum);
     76                void SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum);
    7677                bool InAnalysis(int analysis_type);
    7778                /*}}}*/

  • issm/trunk-jpl/src/c/classes/Loads/Pengrid.cpp

    r16035 r16042  
    204204
    205205        sidlist[0]=node->Sid();
     206}
     207/*}}}*/
     208/*FUNCTION Pengrid::GetNodesLidList{{{*/
     209void Pengrid::GetNodesLidList(int* lidlist){
     210
     211        _assert_(lidlist);
     212        _assert_(node);
     213
     214        lidlist[0]=node->Lid();
    206215}
    207216/*}}}*/
    … …  
    292301/*}}}*/
    293302/*FUNCTION Pengrid::SetwiseNodeConnectivity{{{*/
    294 void Pengrid::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int set1_enum,int set2_enum){
     303void Pengrid::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum){
    295304
    296305        /*Output */
    … …  
    302311                /*flag current node so that no other element processes it*/
    303312                flags[this->node->Lid()]=true;
     313
     314                int counter=0;
     315                while(flagsindices[counter]>=0) counter++;
     316                flagsindices[counter]=this->node->Lid();
    304317
    305318                /*if node is clone, we have an off-diagonal non-zero, else it is a diagonal non-zero*/

  • issm/trunk-jpl/src/c/classes/Loads/Pengrid.h

    r15740 r16042  
    7878                void  CreateJacobianMatrix(Matrix<IssmDouble>* Jff){_error_("Not implemented yet");};
    7979                void  GetNodesSidList(int* sidlist);
     80                void  GetNodesLidList(int* lidlist);
    8081                int   GetNumberOfNodes(void);
    8182                bool  IsPenalty(void);
    … …  
    8384                void  PenaltyCreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* kfs, IssmDouble kmax);
    8485                void  PenaltyCreatePVector(Vector<IssmDouble>* pf, IssmDouble kmax);
    85                 void  SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int set1_enum,int set2_enum);
     86                void  SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum);
    8687                bool  InAnalysis(int analysis_type);
    8788                /*}}}*/

  • issm/trunk-jpl/src/c/classes/Loads/Penpair.cpp

    r16035 r16042  
    151151
    152152        for(int i=0;i<NUMVERTICES;i++) sidlist[i]=nodes[i]->Sid();
     153}
     154/*}}}*/
     155/*FUNCTION Penpair::GetNodesLidList{{{*/
     156void Penpair::GetNodesLidList(int* lidlist){
     157
     158        _assert_(lidlist);
     159        _assert_(nodes);
     160
     161        for(int i=0;i<NUMVERTICES;i++) lidlist[i]=nodes[i]->Lid();
    153162}
    154163/*}}}*/
    … …  
    208217/*}}}*/
    209218/*FUNCTION Penpair::SetwiseNodeConnectivity{{{*/
    210 void Penpair::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int set1_enum,int set2_enum){
     219void Penpair::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum){
    211220
    212221        /*Output */
    … …  
    221230                        /*flag current node so that no other element processes it*/
    222231                        flags[this->nodes[i]->Lid()]=true;
     232
     233                        int counter=0;
     234                        while(flagsindices[counter]>=0) counter++;
     235                        flagsindices[counter]=this->nodes[i]->Lid();
    223236
    224237                        /*if node is clone, we have an off-diagonal non-zero, else it is a diagonal non-zero*/

  • issm/trunk-jpl/src/c/classes/Loads/Penpair.h

    r15771 r16042  
    5858                void  CreateJacobianMatrix(Matrix<IssmDouble>* Jff);
    5959                void  GetNodesSidList(int* sidlist);
     60                void  GetNodesLidList(int* lidlist);
    6061                int   GetNumberOfNodes(void);
    6162                bool  IsPenalty(void);
    … …  
    6364                void  PenaltyCreatePVector(Vector<IssmDouble>* pf, IssmDouble kmax);
    6465                void  PenaltyCreateJacobianMatrix(Matrix<IssmDouble>* Jff,IssmDouble kmax);
    65                 void  SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int set1_enum,int set2_enum);
     66                void  SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum);
    6667                bool  InAnalysis(int analysis_type);
    6768                /*}}}*/

  • issm/trunk-jpl/src/c/classes/Loads/Riftfront.cpp

    r16035 r16042  
    372372}
    373373/*}}}*/
     374/*FUNCTION Riftfront::GetNodesLidList{{{*/
     375void Riftfront::GetNodesLidList(int* lidlist){
     376
     377        _assert_(lidlist);
     378        _assert_(nodes);
     379
     380        for(int i=0;i<NUMVERTICES;i++) lidlist[i]=nodes[i]->Lid();
     381}
     382/*}}}*/
    374383/*FUNCTION Riftfront::GetNumberOfNodes{{{*/
    375384int Riftfront::GetNumberOfNodes(void){
    … …  
    385394/*}}}*/
    386395/*FUNCTION Riftfront::SetwiseNodeConnectivity{{{*/
    387 void Riftfront::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int set1_enum,int set2_enum){
     396void Riftfront::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum){
    388397
    389398        /*Output */
    … …  
    398407                        /*flag current node so that no other element processes it*/
    399408                        flags[this->nodes[i]->Lid()]=true;
     409
     410                        int counter=0;
     411                        while(flagsindices[counter]>=0) counter++;
     412                        flagsindices[counter]=this->nodes[i]->Lid();
    400413
    401414                        /*if node is clone, we have an off-diagonal non-zero, else it is a diagonal non-zero*/

  • issm/trunk-jpl/src/c/classes/Loads/Riftfront.h

    r15771 r16042  
    7979                void  CreateJacobianMatrix(Matrix<IssmDouble>* Jff){_error_("Not implemented yet");};
    8080                void  GetNodesSidList(int* sidlist);
     81                void  GetNodesLidList(int* lidlist);
    8182                int   GetNumberOfNodes(void);
    8283                bool  IsPenalty(void);
    … …  
    8485                void  PenaltyCreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* kfs, IssmDouble kmax);
    8586                void  PenaltyCreatePVector(Vector<IssmDouble>* pf, IssmDouble kmax);
    86                 void  SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int set1_enum,int set2_enum);
     87                void  SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum);
    8788                bool  InAnalysis(int analysis_type);
    8889                /*}}}*/
Note: See TracChangeset for help on using the changeset viewer.