Changeset 4699

Timestamp:

07/21/10 14:47:33 (15 years ago)

Author:

Mathieu Morlighem

Message:

oops... reverted to old version

File:

• issm/trunk/src/c/modules/PenaltyConstraintsx/RiftConstraints.cpp (modified) (1 diff)

issm/trunk/src/c/modules/PenaltyConstraintsx/RiftConstraints.cpp

@@ -1 @@
-s:
+/*!\file RiftConstraints.cpp
+ * \brief: manage penalties for rifts 
+ */
+
+#include "./PenaltyConstraintsLocal.h"
+#include "../../EnumDefinitions/EnumDefinitions.h"
+#include "../../include/include.h"
+#include "../../shared/shared.h"
+
+#define _ZIGZAGCOUNTER_
+
+int RiftConstraints(int* pconverged, int* pnum_unstable_constraints,Loads* loads,int min_mechanical_constraints,int analysis_type){
+
+        int num_unstable_constraints=0;
+        int converged=0;
+        int potential;
+        ISSMERROR(" check analysis_type across all routine!");
+
+        Constrain(&num_unstable_constraints,loads);
+        if(num_unstable_constraints==0)converged=1;
+        
+        
+        if(IsFrozen(loads)){
+                converged=1;
+                num_unstable_constraints=0;
+        }
+        else if(num_unstable_constraints<=min_mechanical_constraints){
+                _printf_("   freezing constraints\n");
+                FreezeConstraints(loads);
+        }
+
+        /*Assign output pointers: */
+        *pconverged=converged;
+        *pnum_unstable_constraints=num_unstable_constraints;
+}
+
+int IsMaterialStable(Loads* loads){
+
+        int i;
+        
+        Riftfront* riftfront=NULL;
+        int found=0;
+        int mpi_found=0;
+
+        /*go though loads, and if non-linearity of the material has converged, let all penalties know: */
+        for (i=0;i<loads->Size();i++){
+
+                if(RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+
+                        if (riftfront->IsMaterialStable()){
+                                found=1;
+                                /*do not break! all penalties should get informed the non-linearity converged!*/
+                        }
+                }
+        }
+
+        #ifdef _PARALLEL_
+        MPI_Reduce (&found,&mpi_found,1,MPI_INT,MPI_SUM,0,MPI_COMM_WORLD );
+        MPI_Bcast(&mpi_found,1,MPI_INT,0,MPI_COMM_WORLD);                
+        found=mpi_found;
+        #endif
+
+        return found;
+}
+
+int RiftIsPresent(Loads* loads,int analysis_type){
+
+
+        int i;
+        
+        int found=0;
+        int mpi_found=0;
+
+        /*go though loads, and figure out if one of the loads is a Riftfront: */
+        for (i=0;i<loads->Size();i++){
+                Load* load=(Load*)loads->GetObjectByOffset(i);
+                if(load->InAnalysis(analysis_type)){
+                        if(RiftfrontEnum==loads->GetEnum(i)){
+                                found=1;
+                                break;
+                        }
+                }
+        }
+
+        #ifdef _PARALLEL_
+        MPI_Reduce (&found,&mpi_found,1,MPI_INT,MPI_SUM,0,MPI_COMM_WORLD );
+        MPI_Bcast(&mpi_found,1,MPI_INT,0,MPI_COMM_WORLD);                
+        found=mpi_found;
+        #endif
+
+        return found;
+}
+
+int IsPreStable(Loads* loads){
+
+
+        int i;
+        
+        Riftfront* riftfront=NULL;
+        int found=0;
+        int mpi_found=0;
+
+        /*go though loads, and figure out if one of the penpair loads is still not stable: */
+        for (i=0;i<loads->Size();i++){
+
+                if(RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+
+                        if (riftfront->PreStable()==0){
+                                found=1;
+                                break;
+                        }
+                }
+        }
+
+        #ifdef _PARALLEL_
+        MPI_Reduce (&found,&mpi_found,1,MPI_INT,MPI_SUM,0,MPI_COMM_WORLD );
+        MPI_Bcast(&mpi_found,1,MPI_INT,0,MPI_COMM_WORLD);                
+        found=mpi_found;
+        #endif
+
+        if (found){
+                /*We found an unstable constraint. : */
+                return 0;
+        }
+        else{
+                return 1;
+        }
+}
+
+int SetPreStable(Loads* loads){
+
+
+        int i;
+        
+        Riftfront* riftfront=NULL;
+        int found=0;
+        int mpi_found=0;
+
+        /*go though loads, and set loads to pre stable.:*/
+        for (i=0;i<loads->Size();i++){
+
+                if(RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+                        riftfront->SetPreStable();
+                }
+        }
+}
+
+int PreConstrain(int* pnum_unstable_constraints,Loads* loads){
+
+        int                     i;
+        
+        /* generic object pointer: */
+        Riftfront* riftfront=NULL;
+
+        int unstable;
+        int sum_num_unstable_constraints;
+        int num_unstable_constraints=0; 
+                
+        /*Enforce constraints: */
+        for (i=0;i<loads->Size();i++){
+
+                if (RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+
+                        riftfront->PreConstrain(&unstable);
+
+                        num_unstable_constraints+=unstable;
+                }
+        }
+
+        #ifdef _PARALLEL_
+        MPI_Reduce (&num_unstable_constraints,&sum_num_unstable_constraints,1,MPI_INT,MPI_SUM,0,MPI_COMM_WORLD );
+        MPI_Bcast(&sum_num_unstable_constraints,1,MPI_INT,0,MPI_COMM_WORLD);                
+        num_unstable_constraints=sum_num_unstable_constraints;
+        #endif
+        
+        /*Assign output pointers: */
+        *pnum_unstable_constraints=num_unstable_constraints;
+
+}
+
+int Constrain(int* pnum_unstable_constraints,Loads* loads){
+
+        int                     i;
+        
+        /* generic object pointer: */
+        Riftfront* riftfront=NULL;
+
+        int unstable;
+        int sum_num_unstable_constraints;
+        int num_unstable_constraints=0; 
+                
+        /*Enforce constraints: */
+        for (i=0;i<loads->Size();i++){
+
+                if (RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+
+                        riftfront->Constrain(&unstable);
+
+                        num_unstable_constraints+=unstable;
+                }
+        }
+
+        #ifdef _PARALLEL_
+        MPI_Reduce (&num_unstable_constraints,&sum_num_unstable_constraints,1,MPI_INT,MPI_SUM,0,MPI_COMM_WORLD );
+        MPI_Bcast(&sum_num_unstable_constraints,1,MPI_INT,0,MPI_COMM_WORLD);                
+        num_unstable_constraints=sum_num_unstable_constraints;
+        #endif
+        
+        /*Assign output pointers: */
+        *pnum_unstable_constraints=num_unstable_constraints;
+
+}
+
+void FreezeConstraints(Loads* loads){
+
+        int                     i;
+        
+        /* generic object pointer: */
+        Riftfront* riftfront=NULL;
+
+        /*Enforce constraints: */
+        for (i=0;i<loads->Size();i++){
+
+                if (RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+
+                        riftfront->FreezeConstraints();
+
+                }
+        }
+
+}
+
+int IsFrozen(Loads* loads){
+
+        int                     i;
+        
+        /* generic object pointer: */
+        Riftfront* riftfront=NULL;
+        int found=0;
+        int mpi_found=0;
+
+        /*Enforce constraints: */
+        for (i=0;i<loads->Size();i++){
+
+                if (RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+                        if (riftfront->IsFrozen()){
+                                found=1;
+                                break;
+                        }
+                }
+        }
+        
+        /*Is there just one found? that would mean we have frozen! : */
+        #ifdef _PARALLEL_
+        MPI_Reduce (&found,&mpi_found,1,MPI_DOUBLE,MPI_MAX,0,MPI_COMM_WORLD );
+        MPI_Bcast(&mpi_found,1,MPI_DOUBLE,0,MPI_COMM_WORLD);                
+        found=mpi_found;
+        #endif
+
+        return found;
+}
+
+int MaxPenetrationInInputs(Loads* loads){
+
+        int                     i;
+        
+        /* generic object pointer: */
+        Riftfront* riftfront=NULL;
+
+        /*rift penetration: */
+        double max_penetration=0;
+        double mpi_max_penetration;
+        double penetration;
+
+        /*Ok, we are going to find the grid pairs which are not penetrating, even though they 
+         * are penalised. We will release only the one with has least <0 penetration. : */
+
+        max_penetration=0;
+        for (i=0;i<loads->Size();i++){
+
+                if (RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+
+                        riftfront->MaxPenetration(&penetration);
+
+                        if (penetration>max_penetration)max_penetration=penetration;
+                }
+        }
+
+        #ifdef _PARALLEL_
+        MPI_Reduce (&max_penetration,&mpi_max_penetration,1,MPI_DOUBLE,MPI_MAX,0,MPI_COMM_WORLD );
+        MPI_Bcast(&mpi_max_penetration,1,MPI_DOUBLE,0,MPI_COMM_WORLD);                
+        max_penetration=mpi_max_penetration;
+        #endif
+
+        /*feed max_penetration to inputs: */
+        for(i=0;i<loads->Size();i++){
+                Load* load=(Load*)loads->GetObjectByOffset(i);
+                load->InputUpdateFromVector(&max_penetration,MaxPenetrationEnum,ConstantEnum);
+        }
+}
+
+int PotentialUnstableConstraints(Loads* loads){
+
+        int                     i;
+        
+        /* generic object pointer: */
+        Riftfront* riftfront=NULL;
+
+        /*Ok, we are going to find the grid pairs which are not penetrating, even though they 
+         * are penalised. We will release only the one with has least <0 penetration. : */
+        int unstable=0;
+        int sum_num_unstable_constraints=0;
+        int num_unstable_constraints=0;
+
+        for (i=0;i<loads->Size();i++){
+
+                if (RiftfrontEnum==loads->GetEnum(i)){
+
+                        riftfront=(Riftfront*)loads->GetObjectByOffset(i);
+
+                        riftfront->PotentialUnstableConstraint(&unstable);
+
+                        num_unstable_constraints+=unstable;
+                }
+        }
+
+        #ifdef _PARALLEL_
+        MPI_Reduce (&num_unstable_constraints,&sum_num_unstable_constraints,1,MPI_INT,MPI_SUM,0,MPI_COMM_WORLD );
+        MPI_Bcast(&sum_num_unstable_constraints,1,MPI_INT,0,MPI_COMM_WORLD);                
+        num_unstable_constraints=sum_num_unstable_constraints;
+        #endif
+
+        return num_unstable_constraints;
+}