source: issm/trunk/src/c/modules/ConstraintsStatex/ConstraintsStatex.cpp@ 16560

/*!\file ConstraintsStatex
 * \brief: set up penalty constraints on loads
 */

#include "./ConstraintsStatex.h"
#include "./ConstraintsStateLocal.h"
#include "../../shared/shared.h"
#include "../../toolkits/toolkits.h"

void ConstraintsStatex(int* pconverged, int* pnum_unstable_constraints,FemModel* femmodel){

        /*output: */
        int converged                     = 1;
        int num_unstable_constraints      = 0;
        int min_mechanical_constraints    = 0;
        int  unstable                     = 0;
        int  sum_num_unstable_constraints = 0;
        int analysis_type;

        /*Display message*/
        if(VerboseModule()) _printf0_("   Constraining penalties\n");

        /*recover parameters: */
        femmodel->parameters->FindParam(&min_mechanical_constraints,StressbalanceRiftPenaltyThresholdEnum);
        femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum);

        /*Rift penalties first*/
        #ifdef _HAVE_RIFTS_
        if(RiftIsPresent(femmodel->loads,analysis_type)){
                RiftConstraintsState(&converged,&num_unstable_constraints,femmodel->loads,min_mechanical_constraints,analysis_type);
        }
        #endif

        /*Deal with pengrid*/
        for(int i=0;i<femmodel->loads->Size();i++){
                Load* load=(Load*)femmodel->loads->GetObjectByOffset(i);
                if(load->InAnalysis(analysis_type)){
                        if(load->ObjectEnum()==PengridEnum){
                                Pengrid* pengrid=(Pengrid*)load;
                                pengrid->ConstraintActivate(&unstable);
                                num_unstable_constraints += unstable;
                        }
                }
        }
        ISSM_MPI_Reduce(&num_unstable_constraints,&sum_num_unstable_constraints,1,ISSM_MPI_INT,ISSM_MPI_SUM,0,IssmComm::GetComm() );
        ISSM_MPI_Bcast(&sum_num_unstable_constraints,1,ISSM_MPI_INT,0,IssmComm::GetComm());                
        num_unstable_constraints=sum_num_unstable_constraints;

        /*Have we converged? : */
        if(num_unstable_constraints) converged=0;

        /*Assign output pointers: */
        *pconverged                = converged;
        *pnum_unstable_constraints = num_unstable_constraints;
}