1 // =============================================================================
2 // PROJECT CHRONO - http://projectchrono.org
3 //
4 // Copyright (c) 2016 projectchrono.org
5 // All right reserved.
6 //
7 // Use of this source code is governed by a BSD-style license that can be found
8 // in the LICENSE file at the top level of the distribution and at
9 // http://projectchrono.org/license-chrono.txt.
10 //
11 // =============================================================================
12 // Authors: Nic Olsen
13 // =============================================================================
14 
15 #include "chrono_distributed/physics/ChDomainDistributed.h"
16 #include "chrono_distributed/physics/ChSystemDistributed.h"
17 
18 #include "chrono_multicore/ChDataManager.h"
19 
20 #include "chrono/core/ChVector.h"
21 #include "chrono/physics/ChBody.h"
22 
23 #include <mpi.h>
24 #include <stdlib.h>
25 #include <iostream>
26 #include <memory>
27 
28 using namespace chrono;
29 
ChDomainDistributed(ChSystemDistributed * sys)30 ChDomainDistributed::ChDomainDistributed(ChSystemDistributed* sys) {
31     this->my_sys = sys;
32     split_axis = 0;
33     split = false;
34     axis_set = false;
35 }
36 
~ChDomainDistributed()37 ChDomainDistributed::~ChDomainDistributed() {}
38 
SetSplitAxis(int i)39 void ChDomainDistributed::SetSplitAxis(int i) {
40     assert(!split);
41     if (i == 0 || i == 1 || i == 2) {
42         split_axis = i;
43         axis_set = true;
44     } else {
45         GetLog() << "Invalid axis\n";
46     }
47 }
48 
SetSimDomain(const ChVector<> & lo,const ChVector<> & hi)49 void ChDomainDistributed::SetSimDomain(const ChVector<>& lo, const ChVector<>& hi) {
50     assert(!split);
51 
52     boxlo = lo;
53     boxhi = hi;
54 
55     double len_x = boxhi.x() - boxlo.x();
56     double len_y = boxhi.y() - boxlo.y();
57     double len_z = boxhi.z() - boxlo.z();
58 
59     if (len_x <= 0 || len_y <= 0 || len_z <= 0)
60         my_sys->ErrorAbort("Invalid domain dimensions.");
61 
62     if (!axis_set) {
63         // Index of the longest domain axis 0=x, 1=y, 2=z
64         split_axis = (len_x >= len_y) ? 0 : 1;
65         split_axis = (len_z >= boxhi[split_axis] - boxlo[split_axis]) ? 2 : split_axis;
66     }
67 
68     SplitDomain();
69 }
70 
SplitDomain()71 void ChDomainDistributed::SplitDomain() {
72     // Length of this subdomain along the long axis
73     double sub_len = (boxhi[split_axis] - boxlo[split_axis]) / my_sys->num_ranks;
74 
75     for (int i = 0; i < 3; i++) {
76         if (split_axis == i) {
77             sublo[i] = boxlo[i] + my_sys->my_rank * sub_len;
78             subhi[i] = sublo[i] + sub_len;
79         } else {
80             sublo[i] = boxlo[i];
81             subhi[i] = boxhi[i];
82         }
83     }
84     split = true;
85 }
86 
GetRank(const ChVector<double> & pos) const87 int ChDomainDistributed::GetRank(const ChVector<double>& pos) const {
88     double sub_len = subhi[split_axis] - sublo[split_axis];
89     double pt = pos[split_axis] - boxlo[split_axis];
90 
91     return (int)(pt / sub_len);
92 }
93 
GetRegion(double pos) const94 distributed::COMM_STATUS ChDomainDistributed::GetRegion(double pos) const {
95     int num_ranks = my_sys->num_ranks;
96     if (num_ranks == 1) {
97         return distributed::OWNED;
98     }
99     int my_rank = my_sys->my_rank;
100     double ghost_layer = my_sys->GetGhostLayer();
101     double high = subhi[split_axis];
102     double low = sublo[split_axis];
103 
104     if (my_rank != 0 && my_rank != num_ranks - 1) {
105         if (pos >= low + ghost_layer && pos < high - ghost_layer) {
106             return distributed::OWNED;
107         } else if (pos >= high && pos < high + ghost_layer) {
108             return distributed::GHOST_UP;
109         } else if (pos >= high - ghost_layer && pos < high) {
110             return distributed::SHARED_UP;
111         } else if (pos >= low && pos < low + ghost_layer) {
112             return distributed::SHARED_DOWN;
113         } else if (pos >= low - ghost_layer && pos < low) {
114             return distributed::GHOST_DOWN;
115         } else if (pos >= high + ghost_layer) {
116             return distributed::UNOWNED_UP;
117         } else if (pos < low - ghost_layer) {
118             return distributed::UNOWNED_DOWN;
119         }
120     }
121 
122     else if (my_rank == 0) {
123         if (pos >= low && pos < high - ghost_layer) {
124             return distributed::OWNED;
125         } else if (pos >= high && pos < high + ghost_layer) {
126             return distributed::GHOST_UP;
127         } else if (pos >= high - ghost_layer && pos < high) {
128             return distributed::SHARED_UP;
129         } else if (pos >= high + ghost_layer) {
130             return distributed::UNOWNED_UP;
131         } else if (pos < low) {
132             return distributed::UNOWNED_DOWN;
133         }
134     }
135 
136     else if (my_rank == num_ranks - 1) {
137         if (pos >= low + ghost_layer && pos < high) {
138             return distributed::OWNED;
139         } else if (pos >= low && pos < low + ghost_layer) {
140             return distributed::SHARED_DOWN;
141         } else if (pos >= low - ghost_layer && pos < low) {
142             return distributed::GHOST_DOWN;
143         } else if (pos >= high) {
144             return distributed::UNOWNED_UP;
145         } else if (pos < low - ghost_layer) {
146             return distributed::UNOWNED_DOWN;
147         }
148     }
149 
150     GetLog() << "Error classifying body\n";
151     return distributed::UNDEFINED;
152 }
153 
GetBodyRegion(int index) const154 distributed::COMM_STATUS ChDomainDistributed::GetBodyRegion(int index) const {
155     return GetRegion(my_sys->data_manager->host_data.pos_rigid[index][split_axis]);
156 }
157 
GetBodyRegion(std::shared_ptr<ChBody> body) const158 distributed::COMM_STATUS ChDomainDistributed::GetBodyRegion(std::shared_ptr<ChBody> body) const {
159     return GetRegion(body->GetPos()[split_axis]);
160 }
161 
PrintDomain()162 void ChDomainDistributed::PrintDomain() {
163     GetLog() << "Domain:\n"
164                 "Box:\n"
165                 "\tX: "
166              << boxlo.x() << " to " << boxhi.x()
167              << "\n"
168                 "\tY: "
169              << boxlo.y() << " to " << boxhi.y()
170              << "\n"
171                 "\tZ: "
172              << boxlo.z() << " to " << boxhi.z()
173              << "\n"
174                 "Subdomain: Rank "
175              << my_sys->my_rank
176              << "\n"
177                 "\tX: "
178              << sublo.x() << " to " << subhi.x()
179              << "\n"
180                 "\tY: "
181              << sublo.y() << " to " << subhi.y()
182              << "\n"
183                 "\tZ: "
184              << sublo.z() << " to " << subhi.z() << "\n";
185 }
186