Code-Eli  0.3.6
newton_raphson_constrained_system_method.hpp
Go to the documentation of this file.
1 /*********************************************************************************
2 * Copyright (c) 2013 David D. Marshall <ddmarsha@calpoly.edu>
3 *
4 * All rights reserved. This program and the accompanying materials
5 * are made available under the terms of the Eclipse Public License v1.0
6 * which accompanies this distribution, and is available at
7 * http://www.eclipse.org/legal/epl-v10.html
8 *
9 * Contributors:
10 * David D. Marshall - initial code and implementation
11 ********************************************************************************/
12 
13 #ifndef eli_mutil_nls_newton_raphson_constrained_system_method_hpp
14 #define eli_mutil_nls_newton_raphson_constrained_system_method_hpp
15 
16 #include <limits>
17 #include <algorithm>
18 
19 #include "eli/code_eli.hpp"
20 
22 
23 namespace eli
24 {
25  namespace mutil
26  {
27  namespace nls
28  {
29  template<typename data__, size_t N__, size_t NSOL__=1>
31  {
32  public:
33  typedef data__ data_type;
35  {
40  };
41 
42  public:
44  : newton_raphson_system_method<data_type, N__, NSOL__>()
45  {
46  for (size_t i=0; i<N__; ++i)
47  {
48  xmin[i]=0;
49  xmax[i]=0;
52  }
53  }
54 
56  : newton_raphson_system_method<data_type, N__, NSOL__>(nrm)
57  {
58  for (size_t i=0; i<N__; ++i)
59  {
60  xmin[i]=nrm.xmin[i];
61  xmax[i]=nrm.xmax[i];
62  xmin_cond[i]=nrm.xmin_cond[i];
63  xmax_cond[i]=nrm.xmax_cond[i];
64  }
65  }
66 
68  {
69  }
70 
72  {
73  for (size_t i=0; i<N__; ++i)
74  {
77  }
78  }
79 
80  void set_periodic_condition(size_t i, const data_type &dmin, const data_type &dmax)
81  {
82  if (i<N__)
83  {
84  xmin[i]=dmin;
85  xmax[i]=dmax;
88  }
89  }
90 
92  {
93  for (size_t i=0; i<N__; ++i)
94  {
96  }
97  }
98  void unset_lower_condition(size_t i)
99  {
100  if (i<N__)
101  {
103  }
104  }
105  void set_lower_condition(size_t i, const data_type &d, end_condition_usage ec)
106  {
107  if (i<N__)
108  {
109  if ( (xmin_cond[i]==NRC_PERIODIC) && (ec!=NRC_PERIODIC) )
110  {
112  }
113  if (ec==NRC_PERIODIC)
114  {
116  }
117 
118  xmin[i]=d;
119  xmin_cond[i]=ec;
120  }
121  }
122  void get_lower_condition(data_type &d, end_condition_usage &ec, size_t i)
123  {
124  if (i<N__)
125  {
126  d=xmin[i];
127  ec=xmin_cond[i];
128  }
129  }
130 
132  {
133  for (size_t i=0; i<N__; ++i)
134  {
136  }
137  }
138  void unset_upper_condition(size_t i)
139  {
140  if (i<N__)
141  {
143  }
144  }
145 
146  void set_upper_condition(size_t i, const data_type &d, end_condition_usage ec)
147  {
148  if (i<N__)
149  {
150  if ( (xmax_cond[i]==NRC_PERIODIC) && (ec!=NRC_PERIODIC) )
151  {
153  }
154  if (ec==NRC_PERIODIC)
155  {
157  }
158 
159  xmax[i]=d;
160  xmax_cond[i]=ec;
161  }
162  }
163  void get_upper_condition(data_type &d, end_condition_usage &ec, size_t i)
164  {
165  if (i<N__)
166  {
167  d=xmax[i];
168  ec=xmax_cond[i];
169  }
170  }
171 
172  private:
176  {
177  size_t i;
179 
180  // Limit constrained values to hit constraint, but allow other terms to remain.
181  for (i=0; i<N__; ++i)
182  {
183  data_type xinew(x(i)+dx(i));
184 
185  // check if min threshold is hit
186  switch(xmin_cond[i])
187  {
188  case(NRC_EXCLUSIVE):
189  {
190  if (xinew<xmin[i])
191  {
192  dx_new(i)=xmin[i]-x(i);
193  }
194  break;
195  }
196  case(NRC_INCLUSIVE):
197  {
198  if (xinew<=xmin[i])
199  {
200  dx_new(i)=(xmin[i]-x(i))*(1-std::numeric_limits<data_type>::epsilon());
201  }
202  break;
203  }
204  default:
205  case(NRC_NOT_USED):
206  {
207  break;
208  }
209  }
210 
211  // check if max threshold is hit
212  switch(xmax_cond[i])
213  {
214  case(NRC_EXCLUSIVE):
215  {
216  if (xinew>xmax[i])
217  {
218  dx_new(i)=xmax[i]-x(i);
219  }
220  break;
221  }
222  case(NRC_INCLUSIVE):
223  {
224  if (xinew>=xmax[i])
225  {
226  dx_new(i)=(xmax[i]-x(i))*(1-std::numeric_limits<data_type>::epsilon());
227  }
228  break;
229  }
230  default:
231  case(NRC_NOT_USED):
232  {
233  break;
234  }
235  }
236  }
237 
238  // Limit the new dx for periodic conditions.
239  for (i=0; i<N__; ++i)
240  {
241  if (xmin_cond[i]==NRC_PERIODIC)
242  {
243  data_type xinew(x[i]+dx_new[i]), period(xmax[i]-xmin[i]);
244 
245  assert(xmax[i]>xmin[i]);
246  assert(period>0);
247 
248  if (xinew<xmin[i])
249  {
250  xinew-=period*std::floor((xinew-xmin[i])/period);
251 
252  assert(xinew>=xmin[i]);
253  assert(xinew<=xmax[i]);
254 
255  dx_new[i]=xinew-x[i];
256  }
257  }
258 
259  if (xmax_cond[i]==NRC_PERIODIC)
260  {
261  data_type xinew(x[i]+dx_new[i]), period(xmax[i]-xmin[i]);
262 
263  assert(xmax[i]>xmin[i]);
264  assert(period>0);
265 
266  if (xinew>xmax[i])
267  {
268  xinew-=period*std::ceil((xinew-xmax[i])/period);
269  dx_new[i]=fmod(xinew, period)-x[i];
270  assert(xinew>=xmin[i]);
271  assert(xinew<=xmax[i]);
272 
273  dx_new[i]=xinew-x[i];
274  }
275  }
276  }
277 
278  return dx_new;
279  }
280 
281  private:
282  data_type xmin[N__], xmax[N__];
284  };
285  }
286  }
287 }
288 #endif
newton_raphson_constrained_system_method()
Definition: newton_raphson_constrained_system_method.hpp:43
virtual iterative_system_root_base< data_type, N__, NSOL__ >::solution_matrix calculate_delta_factor(const typename iterative_system_root_base< data_type, N__, NSOL__ >::solution_matrix &x, const typename iterative_system_root_base< data_type, N__, NSOL__ >::solution_matrix &dx) const
Definition: newton_raphson_constrained_system_method.hpp:174
Definition: math.hpp:20
Definition: newton_raphson_constrained_system_method.hpp:39
end_condition_usage xmin_cond[N__]
Definition: newton_raphson_constrained_system_method.hpp:283
Eigen::Matrix< data__, N__, NSOL__ > solution_matrix
Definition: iterative_system_root_base.hpp:28
end_condition_usage
Definition: newton_raphson_constrained_system_method.hpp:34
Definition: newton_raphson_constrained_system_method.hpp:38
Definition: newton_raphson_constrained_system_method.hpp:30
void unset_lower_condition(size_t i)
Definition: newton_raphson_constrained_system_method.hpp:98
Definition: newton_raphson_constrained_system_method.hpp:37
void unset_upper_condition()
Definition: newton_raphson_constrained_system_method.hpp:131
void set_upper_condition(size_t i, const data_type &d, end_condition_usage ec)
Definition: newton_raphson_constrained_system_method.hpp:146
void unset_lower_condition()
Definition: newton_raphson_constrained_system_method.hpp:91
newton_raphson_constrained_system_method(const newton_raphson_constrained_system_method< data_type, N__, NSOL__ > &nrm)
Definition: newton_raphson_constrained_system_method.hpp:55
void get_lower_condition(data_type &d, end_condition_usage &ec, size_t i)
Definition: newton_raphson_constrained_system_method.hpp:122
end_condition_usage xmax_cond[N__]
Definition: newton_raphson_constrained_system_method.hpp:283
void get_upper_condition(data_type &d, end_condition_usage &ec, size_t i)
Definition: newton_raphson_constrained_system_method.hpp:163
~newton_raphson_constrained_system_method()
Definition: newton_raphson_constrained_system_method.hpp:67
data__ data_type
Definition: newton_raphson_constrained_system_method.hpp:33
Definition: newton_raphson_system_method.hpp:27
data_type xmin[N__]
Definition: newton_raphson_constrained_system_method.hpp:282
Definition: newton_raphson_constrained_system_method.hpp:36
void set_lower_condition(size_t i, const data_type &d, end_condition_usage ec)
Definition: newton_raphson_constrained_system_method.hpp:105
void unset_conditions()
Definition: newton_raphson_constrained_system_method.hpp:71
void set_periodic_condition(size_t i, const data_type &dmin, const data_type &dmax)
Definition: newton_raphson_constrained_system_method.hpp:80
void unset_upper_condition(size_t i)
Definition: newton_raphson_constrained_system_method.hpp:138
data_type xmax[N__]
Definition: newton_raphson_constrained_system_method.hpp:282