Ice Sheet System Model  4.18
Code documentation
Functions
KillIcebergsx.cpp File Reference
#include "./KillIcebergsx.h"
#include "../../shared/shared.h"
#include "../../toolkits/toolkits.h"
#include "../InputUpdateFromVectorx/InputUpdateFromVectorx.h"

Go to the source code of this file.

Functions

void KillIcebergsx (FemModel *femmodel)
 

Function Documentation

◆ KillIcebergsx()

void KillIcebergsx ( FemModel femmodel)

Definition at line 11 of file KillIcebergsx.cpp.

11  {
12 
13  /*Intermediaries*/
14  int lid;
15 
16  /*retrieve vertex info and prepare element flag to speed up process*/
17  int nbv_local = femmodel->vertices->Size();
18  IssmDouble *local_mask = xNewZeroInit<IssmDouble>(nbv_local);
19  bool *element_flag = xNewZeroInit<bool>(femmodel->elements->Size());
20 
21  /*Step 1, go through all elements and put 1 in local_mask if the element is grounded*/
22  for(int i=0;i<femmodel->elements->Size();i++){
23  Element* element=xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(i));
24  if(!element->IsIceInElement()){
25  /*Nothing to do, just flag element to speed up the computation*/
26  element_flag[i] = true;
27  }
28  else{
29  if(element->IsGrounded()){
30  int numvertices = element->GetNumberOfVertices();
31  for(int v=0;v<numvertices;v++) local_mask[element->vertices[v]->Lid()] = 1.;
32  }
33  }
34  }
35 
36  /*Now we have 2 loops, one across cpus, and one for each cpus: we are going
37  * to propagate the mask if an element is connected to a positive mask
38  * already. We then communicate to the other partitions. We stop when the
39  * mask stops changing*/
40  bool keepsyncing = true;
41  while(keepsyncing){
42 
43  /*Get local mask from parallel vector*/
45 
46  /*Local iterations on partition*/
47  bool keepgoing = true;
48  int didsomething = 0;
49  int iter = 1;
50  while(keepgoing){
51  //_printf0_(" -- Kill icebergs: local iteration "<<iter<<"\n");
52 
53  keepgoing = false;
54  didsomething = 0;
55  for(int i=0;i<femmodel->elements->Size();i++){
56  Element* element=xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(i));
57 
58  if(!element_flag[i]){
59  int numvertices = element->GetNumberOfVertices();
60  bool found1 = false;
61  for(int j=0;j<numvertices;j++){
62  lid = element->vertices[j]->Lid();
63  if(local_mask[lid]>0.){
64  found1 = true;
65  break;
66  }
67  }
68  if(found1){
69  element_flag[i] = true;
70  for(int j=0;j<numvertices;j++){
71  lid = element->vertices[j]->Lid();
72  if(local_mask[lid]==0.){
73  local_mask[lid]=1.;
74  keepgoing = true;
75  didsomething = 1;
76  }
77  }
78  }
79  }
80  }
81  iter++;
82  }
83 
84  /*Check how many iterations all cpus did*/
85  int iter_max;
88  if(iter_max==2){
89  /*If iter is only 2, nothing else was changed in the while loop above (iter is initialized as 1 and then ++)*/
90  keepsyncing = false;
91  }
92  }
93 
94  /*Cleanup*/
95  xDelete<bool>(element_flag);
96 
97  /*OK, now deactivate iceberg and count the number of deactivated vertices*/
98  for(int i=0;i<femmodel->elements->Size();i++){
99  Element* element=xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(i));
100 
101  if(element->IsIceInElement()){
102  int numvertices = element->GetNumberOfVertices();
103  bool deactivate = false;
104  for(int j=0;j<numvertices;j++){
105  lid = element->vertices[j]->Lid();
106  if(local_mask[lid]==0.){
107  deactivate = true;
108  break;
109  }
110  }
111 
112  if(deactivate){
113  int numvertices = element->GetNumberOfVertices();
114  IssmDouble* values = xNew<IssmDouble>(numvertices);
115  for(int j=0;j<numvertices;j++) values[j] = 1.; /*Anything >0 = no ice*/
116  element->AddInput2(MaskIceLevelsetEnum,values,P1Enum);
117  xDelete<IssmDouble>(values);
118  }
119  }
120  }
121 
122  /*cleanup*/
123  xDelete<IssmDouble>(local_mask);
124 }
DataSet::Size
int Size()
Definition: DataSet.cpp:399
Element::IsGrounded
bool IsGrounded()
Definition: Element.cpp:2011
FemModel::vertices
Vertices * vertices
Definition: FemModel.h:49
IssmDouble
double IssmDouble
Definition: types.h:37
MaskIceLevelsetEnum
@ MaskIceLevelsetEnum
Definition: EnumDefinitions.h:641
IssmComm::GetComm
static ISSM_MPI_Comm GetComm(void)
Definition: IssmComm.cpp:30
ISSM_MPI_MAX
#define ISSM_MPI_MAX
Definition: issmmpi.h:131
Element::vertices
Vertex ** vertices
Definition: Element.h:49
Vertex::Lid
int Lid(void)
Definition: Vertex.cpp:166
Element
Definition: Element.h:41
P1Enum
@ P1Enum
Definition: EnumDefinitions.h:662
ISSM_MPI_INT
#define ISSM_MPI_INT
Definition: issmmpi.h:127
Element::AddInput2
virtual void AddInput2(int input_enum, IssmDouble *values, int interpolation_enum)
Definition: Element.h:216
FemModel::elements
Elements * elements
Definition: FemModel.h:44
ISSM_MPI_Bcast
int ISSM_MPI_Bcast(void *buffer, int count, ISSM_MPI_Datatype datatype, int root, ISSM_MPI_Comm comm)
Definition: issmmpi.cpp:162
FemModel::SyncLocalVectorWithClonesVerticesAdd
void SyncLocalVectorWithClonesVerticesAdd(IssmDouble *local_vector)
Definition: FemModel.cpp:1405
Element::IsIceInElement
bool IsIceInElement()
Definition: Element.cpp:2021
DataSet::GetObjectByOffset
Object * GetObjectByOffset(int offset)
Definition: DataSet.cpp:334
ISSM_MPI_Reduce
int ISSM_MPI_Reduce(void *sendbuf, void *recvbuf, int count, ISSM_MPI_Datatype datatype, ISSM_MPI_Op op, int root, ISSM_MPI_Comm comm)
Definition: issmmpi.cpp:373
Element::GetNumberOfVertices
virtual int GetNumberOfVertices(void)=0
femmodel
FemModel * femmodel
Definition: esmfbinders.cpp:16