queso-0.57.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
1D1DFunction.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // QUESO - a library to support the Quantification of Uncertainty
5 // for Estimation, Simulation and Optimization
6 //
7 // Copyright (C) 2008-2017 The PECOS Development Team
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the Version 2.1 GNU Lesser General
11 // Public License as published by the Free Software Foundation.
12 //
13 // This library is distributed in the hope that it will be useful,
14 // but WITHOUT ANY WARRANTY; without even the implied warranty of
15 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 // Lesser General Public License for more details.
17 //
18 // You should have received a copy of the GNU Lesser General Public
19 // License along with this library; if not, write to the Free Software
20 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
21 // Boston, MA 02110-1301 USA
22 //
23 //-----------------------------------------------------------------------el-
24 
25 #include <queso/1D1DFunction.h>
26 #include <queso/1DQuadrature.h>
27 
28 namespace QUESO {
29 
30 //*****************************************************
31 // Base 1D->1D class
32 //*****************************************************
34  double minDomainValue,
35  double maxDomainValue)
36  :
37  m_minDomainValue(minDomainValue),
38  m_maxDomainValue(maxDomainValue)
39 {
40  queso_require_less_msg(m_minDomainValue, m_maxDomainValue, "min >= max");
41 }
42 
44 {
45 }
46 
47 double
49 {
50  return m_minDomainValue;
51 }
52 
53 double
55 {
56  return m_maxDomainValue;
57 }
58 
59 double
60 Base1D1DFunction::multiplyAndIntegrate(const Base1D1DFunction& func, unsigned int quadratureOrder, double* resultWithMultiplicationByTAsWell) const
61 {
62  double value = 0.;
63 
64  queso_not_implemented();
65 
66  if (resultWithMultiplicationByTAsWell) { // Just to eliminate INTEL compiler warnings
67  func.value((double) quadratureOrder);
68  }
69 
70  return value;
71 }
72 
73 //*****************************************************
74 // Generic 1D->1D class
75 //*****************************************************
77  double minDomainValue,
78  double maxDomainValue,
79  double (*valueRoutinePtr)(double domainValue, const void* routinesDataPtr),
80  double (*derivRoutinePtr)(double domainValue, const void* routinesDataPtr),
81  const void* routinesDataPtr)
82  :
83  Base1D1DFunction(minDomainValue,maxDomainValue),
84  m_valueRoutinePtr (valueRoutinePtr),
85  m_derivRoutinePtr (derivRoutinePtr),
86  m_routinesDataPtr (routinesDataPtr)
87 {
88 }
89 
91 {
92 }
93 
94 double
95 Generic1D1DFunction::value(double domainValue) const
96 {
97  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
98  std::cerr << "In Generic1D1DFunction::value()"
99  << ": requested x (" << domainValue
100  << ") is out of the interval (" << m_minDomainValue
101  << ", " << m_maxDomainValue
102  << ")"
103  << std::endl;
104  }
105 
106  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
107 
108  return (*m_valueRoutinePtr)(domainValue,m_routinesDataPtr);
109 }
110 
111 double
112 Generic1D1DFunction::deriv(double domainValue) const
113 {
114  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
115  std::cerr << "In Generic1D1DFunction::deriv()"
116  << ": requested x (" << domainValue
117  << ") is out of the interval (" << m_minDomainValue
118  << ", " << m_maxDomainValue
119  << ")"
120  << std::endl;
121  }
122 
123  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
124 
125  return (*m_derivRoutinePtr)(domainValue,m_routinesDataPtr);
126 }
127 
128 //*****************************************************
129 // Constant 1D->1D class
130 //*****************************************************
132  double minDomainValue,
133  double maxDomainValue,
134  double constantValue)
135  :
136  Base1D1DFunction(minDomainValue,maxDomainValue),
137  m_constantValue (constantValue)
138 {
139 }
140 
142 {
143 }
144 
145 double
146 Constant1D1DFunction::value(double domainValue) const
147 {
148  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
149  std::cerr << "In Constant1D1DFunction::value()"
150  << ": requested x (" << domainValue
151  << ") is out of the interval (" << m_minDomainValue
152  << ", " << m_maxDomainValue
153  << ")"
154  << std::endl;
155  }
156 
157  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
158 
159  return m_constantValue;
160 }
161 
162 double
163 Constant1D1DFunction::deriv(double domainValue) const
164 {
165  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
166  std::cerr << "In Constant1D1DFunction::deriv()"
167  << ": requested x (" << domainValue
168  << ") is out of the interval (" << m_minDomainValue
169  << ", " << m_maxDomainValue
170  << ")"
171  << std::endl;
172  }
173 
174  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
175 
176  return 0.;
177 }
178 
179 //*****************************************************
180 // Linear 1D->1D class
181 //*****************************************************
183  double minDomainValue,
184  double maxDomainValue,
185  double referenceDomainValue,
186  double referenceImageValue,
187  double rateValue)
188  :
189  Base1D1DFunction(minDomainValue,maxDomainValue),
190  m_referenceDomainValue (referenceDomainValue),
191  m_referenceImageValue (referenceImageValue),
192  m_rateValue (rateValue)
193 {
194 }
195 
197 {
198 }
199 
200 double
201 Linear1D1DFunction::value(double domainValue) const
202 {
203  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
204  std::cerr << "In Linear1D1DFunction::value()"
205  << ": requested x (" << domainValue
206  << ") is out of the interval (" << m_minDomainValue
207  << ", " << m_maxDomainValue
208  << ")"
209  << std::endl;
210  }
211 
212  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
213 
214  double imageValue = m_referenceImageValue + m_rateValue*(domainValue - m_referenceDomainValue);
215 
216  return imageValue;
217 }
218 
219 double
220 Linear1D1DFunction::deriv(double domainValue) const
221 {
222  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
223  std::cerr << "In Linear1D1DFunction::deriv()"
224  << ": requested x (" << domainValue
225  << ") is out of the interval (" << m_minDomainValue
226  << ", " << m_maxDomainValue
227  << ")"
228  << std::endl;
229  }
230 
231  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
232 
233  return m_rateValue;
234 }
235 
236 //*****************************************************
237 // PiecewiseLinear 1D->1D class
238 //*****************************************************
240  double minDomainValue,
241  double maxDomainValue,
242  const std::vector<double>& referenceDomainValues,
243  double referenceImageValue0,
244  const std::vector<double>& rateValues)
245  :
246  Base1D1DFunction(minDomainValue,maxDomainValue),
247  m_numRefValues (referenceDomainValues.size()),
248  m_referenceDomainValues(referenceDomainValues),
249  m_rateValues (rateValues)
250 {
251  queso_require_not_equal_to_msg(m_numRefValues, 0, "num ref values = 0");
252 
253  queso_require_equal_to_msg(m_numRefValues, rateValues.size(), "num rate values is inconsistent");
254 
255  for (unsigned int i = 1; i < m_numRefValues; ++i) { // Yes, from '1'
256  queso_require_greater_msg(m_referenceDomainValues[i], m_referenceDomainValues[i-1], "reference domain values are inconsistent");
257  }
258 
259  m_referenceImageValues.clear();
260  m_referenceImageValues.resize(m_numRefValues,0.);
261  m_referenceImageValues[0] = referenceImageValue0;
262  for (unsigned int i = 1; i < m_numRefValues; ++i) { // Yes, from '1'
264  }
265 
266  if (false) { // For debug only
267  std::cout << "In PiecewiseLinear1D1DFunction::constructor():"
268  << std::endl;
269  for (unsigned int i = 0; i < m_numRefValues; ++i) {
270  std::cout << "i = " << i
271  << ", m_referenceDomainValues[i] = " << m_referenceDomainValues[i]
272  << ", m_referenceImageValues[i] = " << m_referenceImageValues[i]
273  << ", m_rateValues[i] = " << m_rateValues[i]
274  << std::endl;
275  }
276  }
277 }
278 
280 {
281  m_rateValues.clear();
282  m_referenceImageValues.clear();
283  m_referenceDomainValues.clear();
284 }
285 
286 double
287 PiecewiseLinear1D1DFunction::value(double domainValue) const
288 {
289  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
290  std::cerr << "In PiecewiseLinear1D1DFunction::value()"
291  << ": requested x (" << domainValue
292  << ") is out of the interval (" << m_minDomainValue
293  << ", " << m_maxDomainValue
294  << ")"
295  << std::endl;
296  }
297 
298  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
299 
300  unsigned int i = 0;
301  if (m_numRefValues == 1) {
302  // Nothing else to do
303  }
304  else {
305  bool referenceDomainValueFound = false;
306  while (referenceDomainValueFound == false) {
307  if (domainValue < m_referenceDomainValues[i+1]) {
308  referenceDomainValueFound = true;
309  }
310  else {
311  ++i;
312  if (i == (m_numRefValues-1)) {
313  referenceDomainValueFound = true;
314  }
315  }
316  }
317  }
318  double imageValue = m_referenceImageValues[i] + m_rateValues[i]*(domainValue - m_referenceDomainValues[i]);
319  if (false) { // For debug only
320  std::cout << "In PiecewiseLinear1D1DFunction::value()"
321  << ": domainValue = " << domainValue
322  << ", i = " << i
323  << ", m_referenceDomainValues[i] = " << m_referenceDomainValues[i]
324  << ", m_referenceImageValues[i] = " << m_referenceImageValues[i]
325  << ", m_rateValues[i] = " << m_rateValues[i]
326  << ", imageValue = " << imageValue
327  << std::endl;
328  }
329 
330  return imageValue;
331 }
332 
333 double
334 PiecewiseLinear1D1DFunction::deriv(double domainValue) const
335 {
336  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
337  std::cerr << "In PiecewiseLinear1D1DFunction::deriv()"
338  << ": requested x (" << domainValue
339  << ") is out of the interval (" << m_minDomainValue
340  << ", " << m_maxDomainValue
341  << ")"
342  << std::endl;
343  }
344 
345  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
346 
347  unsigned int i = 0;
348  if (m_numRefValues == 1) {
349  // Nothing else to do
350  }
351  else {
352  bool referenceDomainValueFound = false;
353  while (referenceDomainValueFound == false) {
354  if (domainValue < m_referenceDomainValues[i+1]) {
355  referenceDomainValueFound = true;
356  }
357  else {
358  ++i;
359  queso_require_less_equal_msg(i, m_numRefValues, "too big 'i'");
360  }
361  }
362  }
363 
364  return m_rateValues[i];
365 }
366 
367 //*****************************************************
368 // Quadratic 1D->1D class
369 //*****************************************************
371  double minDomainValue,
372  double maxDomainValue,
373  double a,
374  double b,
375  double c)
376  :
377  Base1D1DFunction(minDomainValue,maxDomainValue),
378  m_a (a),
379  m_b (b),
380  m_c (c)
381 {
382 }
383 
385 {
386 }
387 
388 double
389 Quadratic1D1DFunction::value(double domainValue) const
390 {
391  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
392  std::cerr << "In Quadratic1D1DFunction::value()"
393  << ": requested x (" << domainValue
394  << ") is out of the interval (" << m_minDomainValue
395  << ", " << m_maxDomainValue
396  << ")"
397  << std::endl;
398  }
399 
400  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
401 
402  double imageValue = m_a*domainValue*domainValue + m_b*domainValue + m_c;
403 
404  return imageValue;
405 }
406 
407 double
408 Quadratic1D1DFunction::deriv(double domainValue) const
409 {
410  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
411  std::cerr << "In Quadratic1D1DFunction::deriv()"
412  << ": requested x (" << domainValue
413  << ") is out of the interval (" << m_minDomainValue
414  << ", " << m_maxDomainValue
415  << ")"
416  << std::endl;
417  }
418 
419  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
420 
421  return 2.*m_a*domainValue + m_b;
422 }
423 
424 //*****************************************************
425 // Sampled 1D->1D class
426 //*****************************************************
428  :
429  Base1D1DFunction(-INFINITY,INFINITY)
430 {
431 
432  // UQ_UNAVAILABLE_RANK,
433  // "SampledD1DFunction::deriv()",
434  // "invalid constructor");
435 }
436 
438  const std::vector<double>& domainValues,
439  const std::vector<double>& imageValues)
440  :
441  Base1D1DFunction(domainValues[0],domainValues[domainValues.size()-1]),
442  m_domainValues (domainValues.size(),0.),
443  m_imageValues (imageValues.size(), 0.)
444 {
445  unsigned int tmpSize = m_domainValues.size();
446  for (unsigned int i = 0; i < tmpSize; ++i) {
447  m_domainValues[i] = domainValues[i];
448  m_imageValues [i] = imageValues [i];
449  }
450 }
451 
453 {
454 }
455 
456 double
457 Sampled1D1DFunction::value(double domainValue) const
458 {
459  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
460  std::cerr << "In Sampled1D1DFunction::value()"
461  << ": requested x (" << domainValue
462  << ") is out of the interval (" << m_minDomainValue
463  << ", " << m_maxDomainValue
464  << ")"
465  << std::endl;
466  }
467 
468  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
469 
470  double returnValue = 0.;
471 
472  unsigned int tmpSize = m_domainValues.size();
473  //std::cout << "In Sampled1D1DFunction::value()"
474  // << ": domainValue = " << domainValue
475  // << ", tmpSize = " << tmpSize
476  // << ", m_domainValues[0] = " << m_domainValues[0]
477  // << ", m_domainValues[max] = " << m_domainValues[tmpSize-1]
478  // << std::endl;
479 
480  queso_require_not_equal_to_msg(tmpSize, 0, "m_domainValues.size() = 0");
481 
482  queso_require_greater_equal_msg(domainValue, m_domainValues[0], "domainValue < m_domainValues[0]");
483 
484  queso_require_greater_equal_msg(m_domainValues[tmpSize-1], domainValue, "m_domainValues[max] < domainValue");
485 
486  unsigned int i = 0;
487  for (i = 0; i < tmpSize; ++i) {
488  if (domainValue <= m_domainValues[i]) break;
489  }
490 
491  if (domainValue == m_domainValues[i]) {
492  //if (domainValueWasMatchedExactly) *domainValueWasMatchedExactly = true;
493  returnValue = m_imageValues[i];
494  }
495  else {
496  //if (domainValueWasMatchedExactly) *domainValueWasMatchedExactly = false;
497  double ratio = (domainValue - m_domainValues[i-1])/(m_domainValues[i]-m_domainValues[i-1]);
498  returnValue = m_imageValues[i-1] + ratio * (m_imageValues[i]-m_imageValues[i-1]);
499  }
500 
501  return returnValue;
502 }
503 
504 double
505 Sampled1D1DFunction::deriv(double domainValue) const
506 {
507  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
508  std::cerr << "In Sampled1D1DFunction::deriv()"
509  << ": requested x (" << domainValue
510  << ") is out of the interval (" << m_minDomainValue
511  << ", " << m_maxDomainValue
512  << ")"
513  << std::endl;
514  }
515 
516  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
517 
518  queso_error_msg("this function makes no sense for this class");
519  return 0.;
520 }
521 
522 const std::vector<double>&
524 {
525  return m_domainValues;
526 }
527 
528 const std::vector<double>&
530 {
531  return m_imageValues;
532 }
533 
534 bool
536 {
537  bool result = false;
538 
539  unsigned int tmpSize = m_domainValues.size();
540  for (unsigned int i = 0; i < tmpSize; ++i) {
541  if (domainValue <= m_domainValues[i]) {
542  result = (domainValue == m_domainValues[i]);
543  break;
544  }
545  }
546 
547  return result;
548 }
549 
550 void
552  const std::vector<double>& domainValues,
553  const std::vector<double>& imageValues)
554 {
555  m_domainValues.clear();
556  m_imageValues.clear();
557 
558  unsigned int tmpSize = domainValues.size();
559  m_minDomainValue = domainValues[0];
560  m_maxDomainValue = domainValues[tmpSize-1];
561 
562  m_domainValues.resize(tmpSize,0.);
563  m_imageValues.resize(tmpSize,0.);
564  for (unsigned int i = 0; i < tmpSize; ++i) {
565  m_domainValues[i] = domainValues[i];
566  m_imageValues [i] = imageValues [i];
567  }
568 
569  return;
570 }
571 
572 void
574  const BaseEnvironment& env,
575  std::ofstream& ofsvar,
576  const std::string& prefixName) const
577 {
578  unsigned int tmpSize = m_domainValues.size();
579  if (tmpSize == 0) {
580  tmpSize = 1;
581  ofsvar << "\n" << prefixName << "Time_sub" << env.subIdString() << " = zeros(" << tmpSize << ",1);"
582  << "\n" << prefixName << "Value_sub" << env.subIdString() << " = zeros(" << tmpSize << ",1);";
583  }
584  else {
585  ofsvar << "\n" << prefixName << "Time_sub" << env.subIdString() << " = zeros(" << tmpSize << ",1);"
586  << "\n" << prefixName << "Value_sub" << env.subIdString() << " = zeros(" << tmpSize << ",1);";
587  for (unsigned int i = 0; i < tmpSize; ++i) {
588  ofsvar << "\n" << prefixName << "Time_sub" << env.subIdString() << "(" << i+1 << ",1) = " << m_domainValues[i] << ";"
589  << "\n" << prefixName << "Value_sub" << env.subIdString() << "(" << i+1 << ",1) = " << m_imageValues[i] << ";";
590  }
591  }
592 
593  return;
594 }
595 
596 //*****************************************************
597 // 'ScalarTimesFunc' 1D->1D class
598 //*****************************************************
600  double scalar,
601  const Base1D1DFunction& func)
602  :
603  Base1D1DFunction(func.minDomainValue(),func.maxDomainValue()),
604  m_scalar(scalar),
605  m_func (func)
606 {
607 }
608 
610 {
611 }
612 
613 double
614 ScalarTimesFunc1D1DFunction::value(double domainValue) const
615 {
616  double value = 0.;
617 
618  value += m_scalar*m_func.value(domainValue);
619 
620  return value;
621 }
622 
623 double
624 ScalarTimesFunc1D1DFunction::deriv(double domainValue) const
625 {
626  double value = 0.;
627 
628  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
629  std::cerr << "In ScalarTimes1D1DFunction::deriv()"
630  << ": requested x (" << domainValue
631  << ") is out of the interval (" << m_minDomainValue
632  << ", " << m_maxDomainValue
633  << ")"
634  << std::endl;
635  }
636 
637  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
638 
639  queso_not_implemented();
640 
641  return value;
642 }
643 
644 //*****************************************************
645 // 'FuncTimesFunc' 1D->1D class
646 //*****************************************************
648  const Base1D1DFunction& func1,
649  const Base1D1DFunction& func2)
650  :
651  Base1D1DFunction(std::max(func1.minDomainValue(),func2.minDomainValue()),std::min(func1.maxDomainValue(),func2.maxDomainValue())),
652  m_func1(func1),
653  m_func2(func2)
654 {
655 }
656 
658 {
659 }
660 
661 double
662 FuncTimesFunc1D1DFunction::value(double domainValue) const
663 {
664  double value = 0.;
665 
666  value += m_func1.value(domainValue)*m_func2.value(domainValue);
667 
668  return value;
669 }
670 
671 double
672 FuncTimesFunc1D1DFunction::deriv(double domainValue) const
673 {
674  double value = 0.;
675 
676  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
677  std::cerr << "In FuncTimes1D1DFunction::deriv()"
678  << ": requested x (" << domainValue
679  << ") is out of the interval (" << m_minDomainValue
680  << ", " << m_maxDomainValue
681  << ")"
682  << std::endl;
683  }
684 
685  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
686 
687  queso_not_implemented();
688 
689  return value;
690 }
691 
692 //*****************************************************
693 // 'FuncPlusFunc' 1D->1D class
694 //*****************************************************
696  const Base1D1DFunction& func1,
697  const Base1D1DFunction& func2)
698  :
699  Base1D1DFunction(std::max(func1.minDomainValue(),func2.minDomainValue()),std::min(func1.maxDomainValue(),func2.maxDomainValue())),
700  m_func1(func1),
701  m_func2(func2)
702 {
703 }
704 
706 {
707 }
708 
709 double
710 FuncPlusFunc1D1DFunction::value(double domainValue) const
711 {
712  double value = 0.;
713 
714  value += m_func1.value(domainValue) + m_func2.value(domainValue);
715 
716  return value;
717 }
718 
719 double
720 FuncPlusFunc1D1DFunction::deriv(double domainValue) const
721 {
722  double value = 0.;
723 
724  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
725  std::cerr << "In FuncPlus1D1DFunction::deriv()"
726  << ": requested x (" << domainValue
727  << ") is out of the interval (" << m_minDomainValue
728  << ", " << m_maxDomainValue
729  << ")"
730  << std::endl;
731  }
732 
733  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
734 
735  queso_not_implemented();
736 
737  return value;
738 }
739 
740 //*****************************************************
741 // Lagrange Polynomial 1D->1D class
742 //*****************************************************
744  const std::vector<double>& positionValues,
745  const std::vector<double>* functionValues)
746  :
747  Base1D1DFunction(-INFINITY,INFINITY),
748  m_positionValues(positionValues),
749  m_functionValues(positionValues.size(),1.)
750 {
751  if (functionValues) {
752  queso_require_equal_to_msg(m_positionValues.size(), functionValues->size(), "invalid input");
753  m_functionValues = *functionValues;
754  }
755 }
756 
758 {
759 }
760 
761 double
762 LagrangePolynomial1D1DFunction::value(double domainValue) const
763 {
764  double value = 0.;
765 
766  for (unsigned int k = 0; k < m_positionValues.size(); ++k) {
767  double scaleFactor = 1.;
768  double posK = m_positionValues[k];
769  for (unsigned int j = 0; j < m_positionValues.size(); ++j) {
770  if (j != k) {
771  double posJ = m_positionValues[j];
772  scaleFactor *= (domainValue-posJ)/(posK-posJ);
773  }
774  }
775 
776  //if (m_env.subDisplayFile()) {
777  // *m_env.subDisplayFile() << "In sddTG<K_V,K_M>::lagrange()"
778  // << ": k = " << k
779  // << ", scaleFactor = " << scaleFactor
780  // << std::endl;
781  //}
782 
783  value += scaleFactor * m_functionValues[k];
784 
785  //if (m_env.subDisplayFile()) {
786  // *m_env.subDisplayFile() << "In sddTG<K_V,K_M>::lagrange()"
787  // << ": k = " << k
788  // << ", value = " << value
789  // << std::endl;
790  //}
791  }
792 
793  return value;
794 }
795 
796 double
797 LagrangePolynomial1D1DFunction::deriv(double domainValue) const
798 {
799  double value = 0.;
800 
801  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
802  std::cerr << "In LagrangePolynomial1D1DFunction::deriv()"
803  << ": requested x (" << domainValue
804  << ") is out of the interval (" << m_minDomainValue
805  << ", " << m_maxDomainValue
806  << ")"
807  << std::endl;
808  }
809 
810  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
811 
812  queso_not_implemented();
813 
814  return value;
815 }
816 
817 //*****************************************************
818 // Lagrange Basis 1D->1D class
819 //*****************************************************
821  const std::vector<double>& positionValues,
822  unsigned int basisIndex)
823  :
824  Base1D1DFunction(-INFINITY,INFINITY),
825  m_positionValues(positionValues),
826  m_basisIndex (basisIndex)
827 {
828  queso_require_less_msg(m_basisIndex, m_positionValues.size(), "invalid input");
829 }
830 
832 {
833 }
834 
835 double
836 LagrangeBasis1D1DFunction::value(double domainValue) const
837 {
838  double scaleFactor = 1.;
839 
840  unsigned int k = m_basisIndex;
841  double posK = m_positionValues[k];
842  for (unsigned int j = 0; j < m_positionValues.size(); ++j) {
843  if (j != k) {
844  double posJ = m_positionValues[j];
845  scaleFactor *= (domainValue-posJ)/(posK-posJ);
846  }
847  }
848 
849  return scaleFactor;
850 }
851 
852 double
853 LagrangeBasis1D1DFunction::deriv(double domainValue) const
854 {
855  double value = 0.;
856 
857  if ((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)) {
858  std::cerr << "In LagrangeBasis1D1DFunction::deriv()"
859  << ": requested x (" << domainValue
860  << ") is out of the interval (" << m_minDomainValue
861  << ", " << m_maxDomainValue
862  << ")"
863  << std::endl;
864  }
865 
866  queso_require_msg(!((domainValue < m_minDomainValue) || (domainValue > m_maxDomainValue)), "x out of range");
867 
868  queso_not_implemented();
869 
870  return value;
871 }
872 
873 template <class T>
874 double
876  const ScalarSequence<T>& scalarSeq2,
877  unsigned int initialPos,
878  double scaleValue1,
879  double scaleValue2,
880  const Base1D1DFunction& func1,
881  const Base1D1DFunction& func2,
882  unsigned int quadratureOrder)
883 {
884  double resultValue = 0.;
885 
886  queso_require_equal_to_msg(initialPos, 0, "not implemented yet for initialPos != 0");
887  queso_require_equal_to_msg(scalarSeq1.subSequenceSize(), scalarSeq2.subSequenceSize(), "different sizes");
888 
889  GaussianHermite1DQuadrature quadObj(0.,1.,quadratureOrder);
890  const std::vector<double>& quadPositions = quadObj.positions();
891  const std::vector<double>& quadWeights = quadObj.weights ();
892  queso_require_equal_to_msg(quadPositions.size(), quadWeights.size(), "quadObj has invalid state");
893 
894  unsigned int numQuadraturePositions = quadPositions.size();
895  unsigned int dataSize = scalarSeq1.subSequenceSize();
896  for (unsigned int k = 0; k < dataSize; ++k) {
897  double value1 = 0.;
898  double value2 = 0.;
899  double x1k = scalarSeq1[k];
900  double x2k = scalarSeq2[k];
901  for (unsigned int j = 0; j < numQuadraturePositions; ++j) {
902  value1 += func1.value(scaleValue1*quadPositions[j]+x1k)*quadWeights[j];
903  value2 += func2.value(scaleValue2*quadPositions[j]+x2k)*quadWeights[j];
904  }
905  resultValue += value1*value2;
906  }
907  resultValue *= 1./(2.*M_PI)/((double) dataSize);
908 
909  return resultValue;
910 }
911 
912 } // End namespace QUESO
double deriv(double domainValue) const
Returns the value of the derivative of the linear function at point domainValue.
Definition: 1D1DFunction.C:220
ScalarTimesFunc1D1DFunction(double scalar, const Base1D1DFunction &func)
Default constructor.
Definition: 1D1DFunction.C:599
double deriv(double domainValue) const
TODO: Returns the value of the derivative of the function func1 by another function func2 at the poin...
Definition: 1D1DFunction.C:672
const Base1D1DFunction & m_func2
Definition: 1D1DFunction.h:532
std::vector< double > m_rateValues
Rate value; in , for each =1 .. m_numRefValues.
Definition: 1D1DFunction.h:277
~Constant1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:141
FuncTimesFunc1D1DFunction(const Base1D1DFunction &func1, const Base1D1DFunction &func2)
Default constructor.
Definition: 1D1DFunction.C:647
LagrangePolynomial1D1DFunction(const std::vector< double > &positionValues, const std::vector< double > *functionValues)
Default constructor.
Definition: 1D1DFunction.C:743
double value(double domainValue) const
Returns the value of the constant function at point domainValue.
Definition: 1D1DFunction.C:146
double value(double domainValue) const
Returns the value of the multiplication of a function func1 by another function func2 at the point do...
Definition: 1D1DFunction.C:662
double SubF1F2Gaussian2dKdeIntegral(const ScalarSequence< T > &scalarSeq1, const ScalarSequence< T > &scalarSeq2, unsigned int initialPos, double scaleValue1, double scaleValue2, const Base1D1DFunction &func1, const Base1D1DFunction &func2, unsigned int quadratureOrder)
Calculates the integral of a 2D Gaussian KDE.
Definition: 1D1DFunction.C:875
Quadratic1D1DFunction(double minDomainValue, double maxDomainValue, double a, double b, double c)
Default constructor.
Definition: 1D1DFunction.C:370
double value(double domainValue) const
Returns the value of the (one-dimensional) function at point domainValue, already multiplied by the s...
Definition: 1D1DFunction.C:614
This (virtual) class sets up the environment underlying the use of the QUESO library by an executable...
Definition: Environment.h:198
double value(double domainValue) const
Returns the value of the Lagrange polynomial at point domainValue.
Definition: 1D1DFunction.C:762
virtual double value(double domainValue) const
Returns the value of the sampled function at point domainValue.
Definition: 1D1DFunction.C:457
const void * m_routinesDataPtr
Definition: 1D1DFunction.h:129
virtual void printForMatlab(const BaseEnvironment &env, std::ofstream &ofsvar, const std::string &prefixName) const
Prints the values of the function in Matlab/Octave format.
Definition: 1D1DFunction.C:573
std::vector< double > m_domainValues
Array of the values in the domain of the function (values of the independent variable).
Definition: 1D1DFunction.h:404
virtual double value(double domainValue) const =0
Returns the value of the (one-dimensional) function. See template specialization. ...
std::vector< double > m_imageValues
Array of the values in the image of the function (values of the dependent variable).
Definition: 1D1DFunction.h:407
MonteCarloSGOptions::MonteCarloSGOptions(const BaseEnvironment &env, const char *prefix queso_require_not_equal_to_msg)(m_env.optionsInputFileName(), std::string(""), std::string("this constructor is incompatible with the absence of an options input file"))
PiecewiseLinear1D1DFunction(double minDomainValue, double maxDomainValue, const std::vector< double > &referenceDomainValues, double referenceImageValue0, const std::vector< double > &rateValues)
Default constructor.
Definition: 1D1DFunction.C:239
Sampled1D1DFunction()
Default constructor. It should not be called by the user.
Definition: 1D1DFunction.C:427
~Generic1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:90
std::vector< double > m_positionValues
Definition: 1D1DFunction.h:635
void set(const std::vector< double > &domainValues, const std::vector< double > &imageValues)
Sets the values of the independent (domainValues) and dependent (imageValues) variables of this sampl...
Definition: 1D1DFunction.C:551
const std::vector< double > & domainValues() const
Array of the domain values (values of the independent variable).
Definition: 1D1DFunction.C:523
unsigned int m_numRefValues
Number of points which will be evaluated.
Definition: 1D1DFunction.h:268
double m_b
&#39;Linear&#39; coefficient of the quadratic function; in .
Definition: 1D1DFunction.h:329
double maxDomainValue() const
Returns the maximum value of the domain of the (one-dimensional) function.
Definition: 1D1DFunction.C:54
const std::vector< double > & imageValues() const
Array of the image values (values of the dependent variable).
Definition: 1D1DFunction.C:529
Class for Hermite-Gauss quadrature rule for one-dimensional functions.
Definition: 1DQuadrature.h:194
double deriv(double domainValue) const
TODO: Returns the value of the derivative of the addition of two functions.
Definition: 1D1DFunction.C:720
double value(double domainValue) const
Returns the value of the addition of function func1 and func2 evaluated at the point domainValue...
Definition: 1D1DFunction.C:710
double value(double domainValue) const
Returns the value of the (one-dimensional) function at point domainValue.
Definition: 1D1DFunction.C:95
double minDomainValue() const
Returns the minimum value of the domain of the (one-dimensional) function.
Definition: 1D1DFunction.C:48
double m_referenceImageValue
Reference value in the function image; in .
Definition: 1D1DFunction.h:223
Constant1D1DFunction(double minDomainValue, double maxDomainValue, double constantValue)
Default constructor.
Definition: 1D1DFunction.C:131
Linear1D1DFunction(double minDomainValue, double maxDomainValue, double referenceDomainValue, double referenceImageValue, double rateValue)
Default constructor.
Definition: 1D1DFunction.C:182
double m_referenceDomainValue
Reference value in the function domain; in .
Definition: 1D1DFunction.h:220
Class for handling scalar samples.
double deriv(double domainValue) const
Returns the value of the derivative of the function at point domainValue.
Definition: 1D1DFunction.C:112
double deriv(double domainValue) const
TODO: Returns the value of the derivative of the Lagrange basis at point domainValue.
Definition: 1D1DFunction.C:853
virtual ~Base1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:43
const Base1D1DFunction & m_func2
Definition: 1D1DFunction.h:492
const Base1D1DFunction & m_func
Definition: 1D1DFunction.h:449
Base1D1DFunction(double minDomainValue, double maxDomainValue)
Default constructor.
Definition: 1D1DFunction.C:33
double deriv(double domainValue) const
TODO: Returns the value of the derivative of the function multiplied by the given scalar at point dom...
Definition: 1D1DFunction.C:624
double deriv(double domainValue) const
Returns the value of the derivative of the piecewise-linear function at point domainValue.
Definition: 1D1DFunction.C:334
double(* m_derivRoutinePtr)(double domainValue, const void *routinesDataPtr)
Definition: 1D1DFunction.h:128
std::vector< double > m_functionValues
Definition: 1D1DFunction.h:585
double(* m_valueRoutinePtr)(double domainValue, const void *routinesDataPtr)
Definition: 1D1DFunction.h:127
std::vector< double > m_positionValues
Definition: 1D1DFunction.h:584
double deriv(double domainValue) const
TODO: Returns the value of the derivative of the Lagrange polynomial at point domainValue.
Definition: 1D1DFunction.C:797
double value(double domainValue) const
Returns the value of the quadratic function at point domainValue.
Definition: 1D1DFunction.C:389
LagrangeBasis1D1DFunction(const std::vector< double > &positionValues, unsigned int basisIndex)
Default constructor.
Definition: 1D1DFunction.C:820
std::vector< double > m_referenceDomainValues
Reference values in the piecewise-linear function domain; in , for each =1 .. m_numRefValues.
Definition: 1D1DFunction.h:271
virtual double multiplyAndIntegrate(const Base1D1DFunction &func, unsigned int quadratureOrder, double *resultWithMultiplicationByTAsWell) const
TODO: Multiplies this function with function, and integrates it numerically. See template specializat...
Definition: 1D1DFunction.C:60
const Base1D1DFunction & m_func1
Definition: 1D1DFunction.h:531
double value(double domainValue) const
Returns the value of the linear function at point domainValue.
Definition: 1D1DFunction.C:201
Generic1D1DFunction(double minDomainValue, double maxDomainValue, double(*valueRoutinePtr)(double domainValue, const void *routinesDataPtr), double(*derivRoutinePtr)(double domainValue, const void *routinesDataPtr), const void *routinesDataPtr)
Default constructor.
Definition: 1D1DFunction.C:76
MonteCarloSGOptions::MonteCarloSGOptions(const BaseEnvironment &env, const char *prefix, const McOptionsValues &alternativeOptionsValues queso_require_equal_to_msg)(m_env.optionsInputFileName(), std::string(""), std::string("this constructor is incompatible with the existence of an options input file"))
~Quadratic1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:384
bool domainValueMatchesExactly(double domainValue) const
Checks whether the domain value domainValue matches exactly one of the values in the function domain ...
Definition: 1D1DFunction.C:535
std::vector< double > m_referenceImageValues
Reference values in the piecewise-linear function image; in , for each =1 .. m_numRefValues.
Definition: 1D1DFunction.h:274
const std::string & subIdString() const
Access to the attribute m_subIdString; which stores the string for the sub-environment, and it will be used, for instance, to create the output files for each sub-environment.
Definition: Environment.C:348
double value(double domainValue) const
Returns the value of the piecewise-linear function at point domainValue.
Definition: 1D1DFunction.C:287
double deriv(double domainValue) const
Bogus: Derivative of the function.
Definition: 1D1DFunction.C:505
~FuncTimesFunc1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:657
Class for one-dimensional functions.
Definition: 1D1DFunction.h:53
~Linear1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:196
double value(double domainValue) const
Returns the value of the Lagrange basis at point domainValue.
Definition: 1D1DFunction.C:836
double m_c
&#39;Free&#39; coefficient of the quadratic function; in .
Definition: 1D1DFunction.h:332
virtual ~Sampled1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:452
double deriv(double domainValue) const
Returns the value of the derivative of the quadratic function at point domainValue.
Definition: 1D1DFunction.C:408
double deriv(double domainValue) const
Returns the value of the derivative of the constant function at point domainValue.
Definition: 1D1DFunction.C:163
const Base1D1DFunction & m_func1
Definition: 1D1DFunction.h:491
unsigned int subSequenceSize() const
Size of the sub-sequence of scalars.
double m_a
&#39;Quadratic&#39; coefficient of the quadratic function; in .
Definition: 1D1DFunction.h:326
~LagrangeBasis1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:831
FuncPlusFunc1D1DFunction(const Base1D1DFunction &func1, const Base1D1DFunction &func2)
Default constructor.
Definition: 1D1DFunction.C:695
~FuncPlusFunc1D1DFunction()
Destructor.
Definition: 1D1DFunction.C:705
double m_rateValue
Rate value; in .
Definition: 1D1DFunction.h:226

Generated on Tue Jun 5 2018 19:48:55 for queso-0.57.1 by  doxygen 1.8.5