Mantid
Loading...
Searching...
No Matches
MultiDomainFunction.cpp
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
4// NScD Oak Ridge National Laboratory, European Spallation Source,
5// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6// SPDX - License - Identifier: GPL - 3.0 +
7//----------------------------------------------------------------------
8// Includes
9//----------------------------------------------------------------------
14
15#include <algorithm>
16#include <boost/lexical_cast.hpp>
17#include <set>
18
19namespace Mantid::API {
20
21DECLARE_FUNCTION(MultiDomainFunction)
22
23MultiDomainFunction::MultiDomainFunction() : m_nDomains(0), m_maxIndex(0) { setAttributeValue("NumDeriv", true); }
24
31void MultiDomainFunction::setDomainIndex(size_t funIndex, size_t domainIndex) {
32 m_domains[funIndex] = std::vector<size_t>(1, domainIndex);
34}
35
44void MultiDomainFunction::setDomainIndices(size_t funIndex, const std::vector<size_t> &domainIndices) {
45 m_domains[funIndex] = domainIndices;
47}
48
53 std::set<size_t> dSet;
54 for (auto &domain : m_domains) {
55 if (!domain.second.empty()) {
56 dSet.insert(domain.second.begin(), domain.second.end());
57 }
58 }
59 m_nDomains = dSet.size();
60 m_maxIndex = dSet.empty() ? 0 : *dSet.rbegin();
61}
62
67 m_valueOffsets.clear();
68 m_valueOffsets.emplace_back(0);
69 for (size_t i = 0; i < domain.getNParts(); ++i) {
70 const FunctionDomain &d = domain.getDomain(i);
71 m_valueOffsets.emplace_back(m_valueOffsets.back() + d.size());
72 }
73}
74
77 m_domains.clear();
79}
80
87void MultiDomainFunction::getDomainIndices(size_t funIndex, size_t nDomains, std::vector<size_t> &domains) const {
88 auto it = m_domains.find(funIndex);
89 if (it == m_domains.end()) { // apply to all domains
90 domains.resize(nDomains);
91 for (size_t i = 0; i < domains.size(); ++i) {
92 domains[i] = i;
93 }
94 } else { // apply to selected domains
95 domains.assign(it->second.begin(), it->second.end());
96 }
97}
98
101
106 // works only on CompositeDomain
107 if (!dynamic_cast<const CompositeDomain *>(&domain)) {
108 // make exception if a single member function is defined
109 if (m_maxIndex == 0 && nFunctions() == 1) {
110 getFunction(0)->function(domain, values);
111 return;
112 }
113 throw std::invalid_argument("Non-CompositeDomain passed to MultiDomainFunction.");
114 }
115 const auto &cd = dynamic_cast<const CompositeDomain &>(domain);
116 // domain must not have less parts than m_maxIndex
117 if (cd.getNParts() <= m_maxIndex) {
118 throw std::invalid_argument("CompositeDomain has too few parts (" + std::to_string(cd.getNParts()) +
119 ") for MultiDomainFunction (max index " + std::to_string(m_maxIndex) + ").");
120 }
121 // domain and values must be consistent
122 if (cd.size() != values.size()) {
123 throw std::invalid_argument("MultiDomainFunction: domain and values have different sizes.");
124 }
125
127 // evaluate member functions
128 values.zeroCalculated();
129 for (size_t iFun = 0; iFun < nFunctions(); ++iFun) {
130 // find the domains member function must be applied to
131 std::vector<size_t> domains;
132 getDomainIndices(iFun, cd.getNParts(), domains);
133
134 for (auto &dom : domains) {
135 const FunctionDomain &d = cd.getDomain(dom);
137 getFunction(iFun)->function(d, tmp);
138 values.addToCalculated(m_valueOffsets[dom], tmp);
139 }
140 }
141}
142
145 // works only on CompositeDomain
146 if (!dynamic_cast<const CompositeDomain *>(&domain)) {
147 // make exception if a single member function is defined
148 if (m_maxIndex == 0 && nFunctions() == 1) {
149 getFunction(0)->functionDeriv(domain, jacobian);
150 return;
151 }
152 throw std::invalid_argument("Non-CompositeDomain passed to MultiDomainFunction.");
153 }
154
155 if (getAttribute("NumDeriv").asBool()) {
156 calNumericalDeriv(domain, jacobian);
157 } else {
158 const auto &cd = dynamic_cast<const CompositeDomain &>(domain);
159 // domain must not have less parts than m_maxIndex
160 if (cd.getNParts() < m_maxIndex) {
161 throw std::invalid_argument("CompositeDomain has too few parts (" + std::to_string(cd.getNParts()) +
162 ") for MultiDomainFunction (max index " + std::to_string(m_maxIndex) + ").");
163 }
164
165 jacobian.zero();
167 // evaluate member functions derivatives
168 for (size_t iFun = 0; iFun < nFunctions(); ++iFun) {
169 // find the domains member function must be applied to
170 std::vector<size_t> domains;
171 getDomainIndices(iFun, cd.getNParts(), domains);
172
173 for (auto &dom : domains) {
174 const FunctionDomain &d = cd.getDomain(dom);
175 PartialJacobian J(&jacobian, m_valueOffsets[dom], paramOffset(iFun));
176 getFunction(iFun)->functionDeriv(d, J);
177 }
178 }
179 }
180}
181
187 for (size_t iFun = 0; iFun < nFunctions(); ++iFun) {
188 getFunction(iFun)->iterationStarting();
189 }
190}
191
196 for (size_t iFun = 0; iFun < nFunctions(); ++iFun) {
197 getFunction(iFun)->iterationFinished();
198 }
199}
200
202IFunction::Attribute MultiDomainFunction::getLocalAttribute(size_t funIndex, const std::string &attName) const {
203 if (attName != "domains") {
204 throw std::invalid_argument("MultiDomainFunction does not have attribute " + attName);
205 }
206 if (funIndex >= nFunctions()) {
207 throw std::out_of_range("Function index is out of range.");
208 }
209 try {
210 auto it = m_domains.at(funIndex);
211 if (it.size() == 1 && it.front() == funIndex) {
212 return IFunction::Attribute("i");
213 } else if (!it.empty()) {
214 auto out = std::to_string(it.front());
215 for (auto i = it.begin() + 1; i != it.end(); ++i) {
216 out += "," + std::to_string(*i);
217 }
218 return IFunction::Attribute(out);
219 }
220 } catch (const std::out_of_range &) {
221 return IFunction::Attribute("All");
222 }
223 return IFunction::Attribute("");
224}
225
254void MultiDomainFunction::setLocalAttribute(size_t funIndex, const std::string &attName,
255 const IFunction::Attribute &att) {
256 if (attName != "domains") {
257 throw std::invalid_argument("MultiDomainFunction does not have attribute " + attName);
258 }
259 if (funIndex >= nFunctions()) {
260 throw std::out_of_range("Function index is out of range.");
261 }
262 std::string value = att.asString();
263 auto it = m_domains.find(funIndex);
264
265 if (value == "All") { // fit to all domains
266 if (it != m_domains.end()) {
267 m_domains.erase(it);
268 }
269 return;
270 } else if (value == "i") { // fit to domain with the same index as the function
271 setDomainIndex(funIndex, funIndex);
272 return;
273 } else if (value.empty()) { // do not fit to any domain
274 setDomainIndices(funIndex, std::vector<size_t>());
275 return;
276 }
277
278 // fit to a selection of domains
279 std::vector<size_t> indx;
280 Expression list;
281 list.parse(value);
282 if (list.name() == "+") {
283 if (list.size() != 2 || list.terms()[1].operator_name() != "-") {
284 throw std::runtime_error("MultiDomainFunction: attribute \"domains\" "
285 "expects two integers separated by a \"-\"");
286 }
287 // value looks like "a - b". a and b must be ints and define a range of
288 // domain indices
289 auto start = boost::lexical_cast<size_t>(list.terms()[0].str());
290 size_t end = boost::lexical_cast<size_t>(list.terms()[1].str()) + 1;
291 if (start >= end) {
292 throw std::runtime_error("MultiDomainFunction: attribute \"domains\": wrong range limits.");
293 }
294 indx.resize(end - start);
295 for (size_t i = start; i < end; ++i) {
296 indx[i - start] = i;
297 }
298 } else {
299 // value must be either an int or a list of ints: "a,b,c,..."
300 list.toList();
301 indx.reserve(list.size());
302 std::transform(list.begin(), list.end(), std::back_inserter(indx),
303 [](const Mantid::API::Expression &k) { return boost::lexical_cast<size_t>(k.name()); });
304 }
305 setDomainIndices(funIndex, indx);
306}
307
314std::vector<IFunction_sptr> MultiDomainFunction::createEquivalentFunctions() const {
315 size_t nDomains = m_maxIndex + 1;
316 std::vector<CompositeFunction_sptr> compositeFunctions(nDomains);
317 for (size_t iFun = 0; iFun < nFunctions(); ++iFun) {
318 // find the domains member function must be applied to
319 std::vector<size_t> domains;
320 getDomainIndices(iFun, nDomains, domains);
321
322 for (auto domainIndex : domains) {
323 CompositeFunction_sptr cf = compositeFunctions[domainIndex];
324 if (!cf) {
325 // create a composite function for each domain
327 compositeFunctions[domainIndex] = cf;
328 }
329 // add copies of all functions applied to j-th domain to a single
330 // compositefunction
331 cf->addFunction(FunctionFactory::Instance().createInitialized(getFunction(iFun)->asString()));
332 }
333 }
334 std::vector<IFunction_sptr> outFunctions(nDomains);
335 // fill in the output vector
336 // check functions containing a single member and take it out of the composite
337 for (size_t i = 0; i < compositeFunctions.size(); ++i) {
338 auto fun = compositeFunctions[i];
339 if (!fun || fun->nFunctions() == 0) {
340 throw std::runtime_error("There is no function for domain " + std::to_string(i));
341 }
342 if (fun->nFunctions() > 1) {
343 outFunctions[i] = fun;
344 } else {
345 outFunctions[i] = fun->getFunction(0);
346 }
347 }
348 return outFunctions;
349}
350
351} // namespace Mantid::API
gsl_vector * tmp
double value
The value of the point.
Definition: FitMW.cpp:51
#define DECLARE_FUNCTION(classname)
Macro for declaring a new type of function to be used with the FunctionFactory.
Base class for a composite domain.
virtual const FunctionDomain & getDomain(size_t i) const =0
Return i-th domain.
virtual size_t getNParts() const =0
Return the number of parts in the domain.
size_t paramOffset(size_t i) const
std::size_t nFunctions() const override
Number of functions.
Attribute getAttribute(const std::string &name) const override
Return a value of attribute attName.
CompositeFunction()
Default constructor.
IFunction_sptr getFunction(std::size_t i) const override
Returns the pointer to i-th function.
This class represents an expression made up of names, binary operators and brackets.
Definition: Expression.h:36
std::string name() const
Returns the name of the expression which is a function or variable name.
Definition: Expression.h:71
const std::vector< Expression > & terms() const
Returns the top level terms of the expression (function arguments).
Definition: Expression.h:77
void parse(const std::string &str)
Parse a string and create an expression.
Definition: Expression.cpp:159
iterator end() const
An iterator pointing to the end of the expressions.
Definition: Expression.h:86
iterator begin() const
An iterator pointing to the start of the expressions.
Definition: Expression.h:84
void toList(const std::string &sep=",")
Make sure the expression is a list of expression separated by sep, eg "term1,term2,...
Definition: Expression.cpp:559
size_t size() const
Returns the number of argumens.
Definition: Expression.h:79
Base class that represents the domain of a function.
A class to store values calculated by a function.
void addToCalculated(size_t i, double value)
Add a number to a calculated value.
void zeroCalculated()
Set all calculated values to zero.
size_t size() const
Return the number of values.
Attribute is a non-fitting parameter.
Definition: IFunction.h:282
std::string asString() const
Returns string value if attribute is a string, throws exception otherwise.
Definition: IFunction.cpp:660
std::string asString() const
Writes itself into a string.
Definition: IFunction.cpp:462
void calNumericalDeriv(const FunctionDomain &domain, Jacobian &jacobian)
Calculate numerical derivatives.
Definition: IFunction.cpp:1031
Represents the Jacobian in IFitFunction::functionDeriv.
Definition: Jacobian.h:22
virtual void zero()=0
Zero all matrix elements.
A composite function defined on a CompositeDomain.
size_t getMaxIndex() const
Get the largest domain index.
void setLocalAttribute(size_t funIndex, const std::string &attName, const Attribute &) override
Set a value to attribute attName.
Attribute getLocalAttribute(size_t funIndex, const std::string &attName) const override
Return a value of attribute attName.
void countValueOffsets(const CompositeDomain &domain) const
Count value offsets for each member domain in a CompositeDomain.
void countNumberOfDomains()
Counts number of the domains.
void iterationFinished() override
Called at the end of an iteration.
std::vector< size_t > m_valueOffsets
void functionDeriv(const FunctionDomain &domain, Jacobian &jacobian) override
Derivatives of function with respect to active parameters.
std::map< size_t, std::vector< size_t > > m_domains
Domain index map: finction -> domain.
void getDomainIndices(size_t funIndex, size_t nDomains, std::vector< size_t > &domains) const
Get domain indices for a member function.
size_t getNumberDomains() const override
Get number of domains required by this function.
void iterationStarting() override
Called at the start of each iteration.
void setDomainIndices(size_t funIndex, const std::vector< size_t > &domainIndices)
Associate a function and a list of domains.
void clearDomainIndices()
Clear all domain indices.
void setDomainIndex(size_t funIndex, size_t domainIndex)
Associate a function and a domain.
size_t m_nDomains
Number of domains this MultiDomainFunction operates on.
std::vector< IFunction_sptr > createEquivalentFunctions() const override
Create a list of equivalent functions.
size_t m_maxIndex
Maximum domain index.
void function(const FunctionDomain &domain, FunctionValues &values) const override
Function you want to fit to.
A Jacobian for individual functions.
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
std::shared_ptr< CompositeFunction > CompositeFunction_sptr
shared pointer to the composite function base class
std::string to_string(const wide_integer< Bits, Signed > &n)