/[escript]/branches/diaplayground/ripley/src/Brick.cpp
ViewVC logotype

Diff of /branches/diaplayground/ripley/src/Brick.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 3913 by caltinay, Tue Jun 19 06:21:58 2012 UTC revision 3962 by jfenwick, Tue Sep 11 09:12:07 2012 UTC
# Line 32  namespace ripley { Line 32  namespace ripley {
32  Brick::Brick(int n0, int n1, int n2, double x0, double y0, double z0,  Brick::Brick(int n0, int n1, int n2, double x0, double y0, double z0,
33               double x1, double y1, double z1, int d0, int d1, int d2) :               double x1, double y1, double z1, int d0, int d1, int d2) :
34      RipleyDomain(3),      RipleyDomain(3),
     m_gNE0(n0),  
     m_gNE1(n1),  
     m_gNE2(n2),  
35      m_x0(x0),      m_x0(x0),
36      m_y0(y0),      m_y0(y0),
37      m_z0(z0),      m_z0(z0),
38      m_l0(x1-x0),      m_l0(x1-x0),
39      m_l1(y1-y0),      m_l1(y1-y0),
40      m_l2(z1-z0),      m_l2(z1-z0)
     m_NX(d0),  
     m_NY(d1),  
     m_NZ(d2)  
41  {  {
42        // ignore subdivision parameters for serial run
43        if (m_mpiInfo->size == 1) {
44            d0=1;
45            d1=1;
46            d2=1;
47        }
48    
49        bool warn=false;
50        // if number of subdivisions is non-positive, try to subdivide by the same
51        // ratio as the number of elements
52        if (d0<=0 && d1<=0 && d2<=0) {
53            warn=true;
54            d0=(int)(pow(m_mpiInfo->size*(n0+1)*(n0+1)/((float)(n1+1)*(n2+1)), 1.f/3));
55            d1=(int)(d0*n1/(float)n0);
56            d2=m_mpiInfo->size/(d0*d1);
57            if (d0*d1*d2 != m_mpiInfo->size) {
58                // ratios not the same so leave "smallest" side undivided and try
59                // dividing 2 sides only
60                if (n0>=n1) {
61                    if (n1>=n2) {
62                        d0=d1=0;
63                        d2=1;
64                    } else {
65                        d0=d2=0;
66                        d1=1;
67                    }
68                } else {
69                    if (n0>=n2) {
70                        d0=d1=0;
71                        d2=1;
72                    } else {
73                        d0=1;
74                        d1=d2=0;
75                    }
76                }
77            }
78        }
79        if (d0<=0 && d1<=0) {
80            warn=true;
81            d0=int(sqrt(m_mpiInfo->size*(n0+1)/(float)(n1+1)));
82            d1=m_mpiInfo->size/d0;
83            if (d0*d1*d2 != m_mpiInfo->size) {
84                // ratios not the same so subdivide side with more elements only
85                if (n0>n1) {
86                    d0=0;
87                    d1=1;
88                } else {
89                    d0=1;
90                    d1=0;
91                }
92            }
93        } else if (d0<=0 && d2<=0) {
94            warn=true;
95            d0=int(sqrt(m_mpiInfo->size*(n0+1)/(float)(n2+1)));
96            d2=m_mpiInfo->size/d0;
97            if (d0*d1*d2 != m_mpiInfo->size) {
98                // ratios not the same so subdivide side with more elements only
99                if (n0>n2) {
100                    d0=0;
101                    d2=1;
102                } else {
103                    d0=1;
104                    d2=0;
105                }
106            }
107        } else if (d1<=0 && d2<=0) {
108            warn=true;
109            d1=int(sqrt(m_mpiInfo->size*(n1+1)/(float)(n2+1)));
110            d2=m_mpiInfo->size/d1;
111            if (d0*d1*d2 != m_mpiInfo->size) {
112                // ratios not the same so subdivide side with more elements only
113                if (n1>n2) {
114                    d1=0;
115                    d2=1;
116                } else {
117                    d1=1;
118                    d2=0;
119                }
120            }
121        }
122        if (d0<=0) {
123            // d1,d2 are preset, determine d0
124            d0=m_mpiInfo->size/(d1*d2);
125        } else if (d1<=0) {
126            // d0,d2 are preset, determine d1
127            d1=m_mpiInfo->size/(d0*d2);
128        } else if (d2<=0) {
129            // d0,d1 are preset, determine d2
130            d2=m_mpiInfo->size/(d0*d1);
131        }
132    
133        m_NX=d0;
134        m_NY=d1;
135        m_NZ=d2;
136    
137      // ensure number of subdivisions is valid and nodes can be distributed      // ensure number of subdivisions is valid and nodes can be distributed
138      // among number of ranks      // among number of ranks
139      if (m_NX*m_NY*m_NZ != m_mpiInfo->size)      if (m_NX*m_NY*m_NZ != m_mpiInfo->size)
140          throw RipleyException("Invalid number of spatial subdivisions");          throw RipleyException("Invalid number of spatial subdivisions");
141    
142      if ((n0+1)%m_NX > 0 || (n1+1)%m_NY > 0 || (n2+1)%m_NZ > 0)      if (warn) {
143          throw RipleyException("Number of elements+1 must be separable into number of ranks in each dimension");          cout << "Warning: Automatic domain subdivision (d0=" << d0 << ", d1="
144                << d1 << ", d2=" << d2 << "). This may not be optimal!" << endl;
145        }
146    
147        if ((n0+1)%m_NX > 0) {
148            double Dx=m_l0/n0;
149            n0=(int)round((float)(n0+1)/d0+0.5)*d0-1;
150            m_l0=Dx*n0;
151            cout << "Warning: Adjusted number of elements and length. N0="
152                << n0 << ", l0=" << m_l0 << endl;
153        }
154        if ((n1+1)%m_NY > 0) {
155            double Dy=m_l1/n1;
156            n1=(int)round((float)(n1+1)/d1+0.5)*d1-1;
157            m_l1=Dy*n1;
158            cout << "Warning: Adjusted number of elements and length. N1="
159                << n1 << ", l1=" << m_l1 << endl;
160        }
161        if ((n2+1)%m_NZ > 0) {
162            double Dz=m_l2/n2;
163            n2=(int)round((float)(n2+1)/d2+0.5)*d2-1;
164            m_l2=Dz*n2;
165            cout << "Warning: Adjusted number of elements and length. N2="
166                << n2 << ", l2=" << m_l2 << endl;
167        }
168    
169        m_gNE0=n0;
170        m_gNE1=n1;
171        m_gNE2=n2;
172    
173      if ((m_NX > 1 && (n0+1)/m_NX<2) || (m_NY > 1 && (n1+1)/m_NY<2) || (m_NZ > 1 && (n2+1)/m_NZ<2))      if ((m_NX > 1 && (n0+1)/m_NX<2) || (m_NY > 1 && (n1+1)/m_NY<2) || (m_NZ > 1 && (n2+1)/m_NZ<2))
174          throw RipleyException("Too few elements for the number of ranks");          throw RipleyException("Too few elements for the number of ranks");

Legend:
Removed from v.3913  
changed lines
  Added in v.3962

  ViewVC Help
Powered by ViewVC 1.1.26