Changeset 24630


Ignore:
Timestamp:
03/03/20 12:19:15 (5 years ago)
Author:
Mathieu Morlighem
Message:

CHG: speeding up ConstraintState if no penalties and no rift found

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

Legend:

Unmodified
Added
Removed

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

    r23587 r24630  
    1919#include "../../shared/Enum/EnumDefinitions.h"
    2020#include "../../shared/Exceptions/exceptions.h"
     21#include "../../shared/MemOps/MemOps.h"
     22#include "../../shared/io/Marshalling/Marshalling.h"
    2123#include "./Loads.h"
    2224#include "./Load.h"
    … …  
    2729/*Object constructors and destructor*/
    2830Loads::Loads(){/*{{{*/
    29         enum_type=LoadsEnum;
     31        this->enum_type=LoadsEnum;
     32        this->numrifts     = 0;
     33        this->numpenalties = 0;
    3034        return;
    3135}
    … …  
    3640/*}}}*/
    3741
     42Loads* Loads::Copy() {/*{{{*/
     43
     44        int num_proc = IssmComm::GetSize();
     45
     46        /*Copy dataset*/
     47        Loads* output=new Loads();
     48        output->sorted=this->sorted;
     49        output->numsorted=this->numsorted;
     50        output->presorted=this->presorted;
     51        for(vector<Object*>::iterator obj=this->objects.begin() ; obj < this->objects.end(); obj++ ){
     52                output->AddObject((*obj)->copy());
     53        }
     54
     55        /*Build id_offsets and sorted_ids*/
     56        int objsize = this->numsorted;
     57        output->id_offsets=NULL;
     58        output->sorted_ids=NULL;
     59        if(this->sorted && objsize>0 && this->id_offsets){     
     60                output->id_offsets=xNew<int>(objsize);
     61                xMemCpy<int>(output->id_offsets,this->id_offsets,objsize);
     62        }
     63        if(this->sorted && objsize>0 && this->sorted_ids){
     64                output->sorted_ids=xNew<int>(objsize);
     65                xMemCpy<int>(output->sorted_ids,this->sorted_ids,objsize);
     66        }
     67
     68        /*Copy other fields*/
     69        output->numrifts = this->numrifts;
     70        output->numpenalties = this->numpenalties;
     71
     72        return output;
     73}
     74/*}}}*/
     75void  Loads::Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction){ /*{{{*/
     76
     77        int num_procs=IssmComm::GetSize();
     78        int test = num_procs;
     79        MARSHALLING_ENUM(LoadsEnum);
     80        MARSHALLING(numrifts);
     81        MARSHALLING(numpenalties);
     82
     83        DataSet::Marshall(pmarshalled_data,pmarshalled_data_size,marshall_direction);
     84}
     85/*}}}*/
     86
    3887/*Numerics:*/
    3988void Loads::Configure(Elements* elements,Loads* loads, Nodes* nodes, Vertices* vertices, Materials* materials,Parameters* parameters){/*{{{*/
    4089
    4190        vector<Object*>::iterator object;
    42         Load* load=NULL;
     91        for(object=objects.begin() ; object < objects.end(); object++){
     92                Load* load=xDynamicCast<Load*>(*object);
     93                load->Configure(elements,loads,nodes,vertices,materials,parameters);
     94        }
     95}
     96/*}}}*/
     97void Loads::Finalize(){/*{{{*/
    4398
    44         for ( object=objects.begin() ; object < objects.end(); object++ ){
     99        /*Count Rifts and penalties*/
     100        int ispenalty=0;
     101        int isrift=0;
     102        int allcount;
    45103
    46                 load=xDynamicCast<Load*>(*object);
    47                 load->Configure(elements,loads,nodes,vertices,materials,parameters);
     104        /*Now go through all loads, and get how many nodes they own, unless they are clone nodes: */
     105        for(int i=0;i<this->Size();i++){
     106                Load* load=xDynamicCast<Load*>(this->GetObjectByOffset(i));
     107                if(load->IsPenalty()){
     108                        ispenalty++;
     109                }
     110      if(load->ObjectEnum()==RiftfrontEnum){
     111         isrift++;
     112      }
     113        }
    48114
    49         }
     115        /*Grab sum of all cpus: */
     116        ISSM_MPI_Allreduce((void*)&ispenalty,(void*)&allcount,1,ISSM_MPI_INT,ISSM_MPI_SUM,IssmComm::GetComm());
     117        this->numpenalties = allcount;
     118
     119        ISSM_MPI_Allreduce((void*)&isrift,(void*)&allcount,1,ISSM_MPI_INT,ISSM_MPI_SUM,IssmComm::GetComm());
     120        this->numrifts= allcount;
    50121
    51122}
    … …  
    53124bool Loads::IsPenalty(){/*{{{*/
    54125
    55         int ispenalty=0;
    56         int allispenalty=0;
    57 
    58         /*Now go through all loads, and get how many nodes they own, unless they are clone nodes: */
    59         for(int i=0;i<this->Size();i++){
    60                 Load* load=xDynamicCast<Load*>(this->GetObjectByOffset(i));
    61                 if(load->IsPenalty()) ispenalty++;
    62         }
    63 
    64         /*Grab sum of all cpus: */
    65         ISSM_MPI_Allreduce((void*)&ispenalty,(void*)&allispenalty,1,ISSM_MPI_INT,ISSM_MPI_SUM,IssmComm::GetComm());
    66         ispenalty=allispenalty;
    67 
    68         if(ispenalty)
     126        if(this->numpenalties>0)
    69127         return true;
    70128        else

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

    r23587 r24630  
    1818        public:
    1919
     20                int numrifts;
     21                int numpenalties;
     22
    2023                /*constructors, destructors*/
    2124                Loads();
    2225                ~Loads();
    2326
     27                /*Objects virtual functions*/
     28                Loads* Copy();
     29                void   Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction);
     30
    2431                /*numerics*/
    2532                void  Configure(Elements* elements,Loads* loads, Nodes* nodes, Vertices* vertices, Materials* materials,Parameters* parameters);
    2633                bool  IsPenalty();
     34                void  Finalize();
    2735                int   MaxNumNodes();
    2836                int   NumberOfLoads();

  • TabularUnified issm/trunk-jpl/src/c/modules/ConstraintsStatex/ConstraintsStatex.cpp

    r23588 r24630  
    1010void ConstraintsStatex(int* pconverged, int* pnum_unstable_constraints,FemModel* femmodel){
    1111
     12        /*Early return if no rift and no penalties*/
     13        if(femmodel->loads->numrifts == 0 && femmodel->loads->numpenalties == 0){
     14                *pconverged                = 0;
     15                *pnum_unstable_constraints = 0;
     16                return;
     17        }
     18
    1219        /*output: */
    13         int converged                     = 1;
    14         int num_unstable_constraints      = 0;
    15         int min_mechanical_constraints    = 0;
    16         int  unstable                     = 0;
    17         int  sum_num_unstable_constraints = 0;
     20        int converged                    = 1;
     21        int num_unstable_constraints     = 0;
     22        int min_mechanical_constraints   = 0;
     23        int unstable                     = 0;
     24        int sum_num_unstable_constraints = 0;
    1825        int analysis_type;
    1926
    … …  
    2633
    2734        /*Rift penalties first*/
    28         if(RiftIsPresent(femmodel->loads,analysis_type)){
     35        if(femmodel->loads->numrifts){
    2936                RiftConstraintsState(&converged,&num_unstable_constraints,femmodel->loads,min_mechanical_constraints,analysis_type);
    3037        }
    3138
    3239        /*Deal with pengrid*/
    33         for(int i=0;i<femmodel->loads->Size();i++){
    34                 Load* load=(Load*)femmodel->loads->GetObjectByOffset(i);
    35                 if(load->ObjectEnum()==PengridEnum){
    36                         Pengrid* pengrid=(Pengrid*)load;
    37                         pengrid->ConstraintActivate(&unstable);
    38                         num_unstable_constraints += unstable;
     40        if(femmodel->loads->numpenalties){
     41                for(int i=0;i<femmodel->loads->Size();i++){
     42                        Load* load=(Load*)femmodel->loads->GetObjectByOffset(i);
     43                        if(load->ObjectEnum()==PengridEnum){
     44                                Pengrid* pengrid=(Pengrid*)load;
     45                                pengrid->ConstraintActivate(&unstable);
     46                                num_unstable_constraints += unstable;
     47                        }
    3948                }
    4049        }
     50
    4151        ISSM_MPI_Reduce(&num_unstable_constraints,&sum_num_unstable_constraints,1,ISSM_MPI_INT,ISSM_MPI_SUM,0,IssmComm::GetComm() );
    4252        ISSM_MPI_Bcast(&sum_num_unstable_constraints,1,ISSM_MPI_INT,0,IssmComm::GetComm());               

  • TabularUnified issm/trunk-jpl/src/c/modules/ConstraintsStatex/ConstraintsStatex.h

    r16206 r24630  
    99
    1010/* local prototypes: */
    11 int  RiftIsPresent(Loads* loads,int analysis_type);
    1211void ConstraintsStatex(int* pconverged, int* pnum_unstable_constraints,FemModel* femmodel);
    1312

  • TabularUnified issm/trunk-jpl/src/c/modules/ConstraintsStatex/RiftConstraintsState.cpp

    r24623 r24630  
    66
    77/*current module: */
    8 /*RiftIsPresent(Loads* loads,int configuration_type){{{*/
    9 int RiftIsPresent(Loads* loads,int configuration_type){
    10 
    11         int i;
    12 
    13         int found=0;
    14         int mpi_found=0;
    15 
    16         /*go though loads, and figure out if one of the loads is a Riftfront: */
    17         for (i=0;i<loads->Size();i++){
    18                 Load* load=(Load*)loads->GetObjectByOffset(i);
    19                 if(RiftfrontEnum==loads->GetEnum(i)){
    20                         found=1;
    21                         break;
    22                 }
    23         }
    24 
    25         ISSM_MPI_Reduce (&found,&mpi_found,1,ISSM_MPI_INT,ISSM_MPI_SUM,0,IssmComm::GetComm() );
    26         ISSM_MPI_Bcast(&mpi_found,1,ISSM_MPI_INT,0,IssmComm::GetComm());               
    27         found=mpi_found;
    28 
    29         return found;
    30 }
    31 /*}}}*/
    328/*RiftConstraintsState(int* pconverged, int* pnum_unstable_constraints,Loads* loads,int min_mechanical_constraints,int configuration_type){{{*/
    339void RiftConstraintsState(int* pconverged, int* pnum_unstable_constraints,Loads* loads,int min_mechanical_constraints,int configuration_type){

  • TabularUnified issm/trunk-jpl/src/c/modules/ModelProcessorx/ModelProcessorx.cpp

    r24335 r24630  
    6969                loads[i]->Presort();
    7070                nodes[i]->Presort();
     71
     72                /*Finalize loads (count pengrids,penpairs,rifts,etc)*/
     73                loads[i]->Finalize();
    7174        }
    7275

  • TabularUnified issm/trunk-jpl/src/c/modules/ResetConstraintsx/ResetConstraintsx.cpp

    r23588 r24630  
    2424
    2525        /*Deal with rift first*/
    26         if(RiftIsPresent(femmodel->loads,analysis_type)){
     26        if(femmodel->loads->numrifts){
    2727                _error_("rift constraints reset not supported yet!");
    2828        }
    … …  
    4343
    4444        /*Deal with rift first*/
    45         if(RiftIsPresent(femmodel->loads,femmodel->analysis_type_list[femmodel->analysis_counter])){
     45        if(femmodel->loads->numrifts){
    4646                _error_("rift constraints reset not supported yet!");
    4747        }
Note: See TracChangeset for help on using the changeset viewer.