1 // -*- C++ -*-
2 /**
3 * @brief LogNormal distribution with mu and sigma over mu as parameters
4 *
5 * Copyright 2005-2021 Airbus-EDF-IMACS-ONERA-Phimeca
6 *
7 * This library is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU Lesser General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * This library is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU Lesser General Public License for more details.
16 *
17 * You should have received a copy of the GNU Lesser General Public License
18 * along with this library. If not, see <http://www.gnu.org/licenses/>.
19 *
20 */
21
22 #include "openturns/LogNormalFactory.hxx"
23 #include "openturns/SquareMatrix.hxx"
24 #include "openturns/LogNormalMuSigma.hxx"
25 #include "openturns/LogNormalMuSigmaOverMu.hxx"
26 #include "openturns/PersistentObjectFactory.hxx"
27
28 BEGIN_NAMESPACE_OPENTURNS
29
30 CLASSNAMEINIT(LogNormalMuSigmaOverMu)
31 static const Factory<LogNormalMuSigmaOverMu> Factory_LogNormalMuSigmaOverMu;
32
33 /* Default constructor */
LogNormalMuSigmaOverMu()34 LogNormalMuSigmaOverMu::LogNormalMuSigmaOverMu()
35 : DistributionParametersImplementation()
36 , mu_(exp(0.5))
37 , sigmaOverMu_(sqrt(exp(2.0) - exp(1.0)) / exp(0.5))
38 , gamma_(0.)
39 {
40 // Nothing to do
41 }
42
LogNormalMuSigmaOverMu(const Scalar mu,const Scalar sigmaOverMu,const Scalar gamma)43 LogNormalMuSigmaOverMu::LogNormalMuSigmaOverMu(const Scalar mu, const Scalar sigmaOverMu, const Scalar gamma)
44 : DistributionParametersImplementation()
45 , mu_(mu)
46 , sigmaOverMu_(sigmaOverMu)
47 , gamma_(gamma)
48 {
49 if (mu == 0.0) throw InvalidArgumentException(HERE) << "mu cannot be null in the parameter set (mu, sigmaOverMu)";
50 if (!(sigmaOverMu * mu > 0.0)) throw InvalidArgumentException(HERE) << "sigmaOverMu*mu must be > 0, here sigmaOverMu*mu=" << sigmaOverMu*mu;
51 if (mu <= gamma) throw InvalidArgumentException(HERE) << "mu must be greater than gamma, here mu=" << mu << " and gamma=" << gamma;
52 }
53
54 /* Virtual constructor */
clone() const55 LogNormalMuSigmaOverMu * LogNormalMuSigmaOverMu::clone() const
56 {
57 return new LogNormalMuSigmaOverMu(*this);
58 }
59
60 /* Comparison operator */
operator ==(const LogNormalMuSigmaOverMu & other) const61 Bool LogNormalMuSigmaOverMu::operator ==(const LogNormalMuSigmaOverMu & other) const
62 {
63 return (this == &other);
64 }
65
66
67 /* Build a distribution based on a set of native parameters */
getDistribution() const68 Distribution LogNormalMuSigmaOverMu::getDistribution() const
69 {
70 Point newParameters(3);
71 newParameters[0] = mu_;
72 newParameters[1] = sigmaOverMu_;
73 newParameters[2] = gamma_;
74
75 Point nativeParameters(operator()(newParameters));
76
77 return LogNormalFactory().build(nativeParameters);
78 }
79
80
81 /* Compute jacobian / native parameters */
gradient() const82 Matrix LogNormalMuSigmaOverMu::gradient() const
83 {
84 // compute the jacobian of the transformation (mulog, sigmalog, gamma) -> (mu, sigma, gamma)
85 const LogNormalMuSigma muSigmaParameters(mu_, sigmaOverMu_ * mu_, gamma_);
86 const Matrix muSigmaJacobian = muSigmaParameters.gradient();
87
88 // compute the jacobian of the transformation (mu, sigma, gamma) -> (mu, sigma/mu, gamma)
89 SquareMatrix muSigmaOverMuJacobian(IdentityMatrix(3));
90 muSigmaOverMuJacobian(0, 1) = sigmaOverMu_;
91 muSigmaOverMuJacobian(1, 1) = mu_;
92
93 return muSigmaOverMuJacobian * muSigmaJacobian;
94 }
95
96
97 /* Conversion operator */
operator ()(const Point & inP) const98 Point LogNormalMuSigmaOverMu::operator () (const Point & inP) const
99 {
100 if (inP.getDimension() != 3) throw InvalidArgumentException(HERE) << "the given point must have dimension=3, here dimension=" << inP.getDimension();
101 const Scalar mu = inP[0];
102 const Scalar sigmaOverMu = inP[1];
103 const Scalar gamma = inP[2];
104
105 if (mu == 0.0) throw InvalidArgumentException(HERE) << "mu cannot be null in the parameter set (mu, sigmaOverMu)";
106 if (!(sigmaOverMu * mu > 0.0)) throw InvalidArgumentException(HERE) << "sigmaOverMu*mu must be > 0, here sigmaOverMu*mu=" << sigmaOverMu*mu;
107 if (mu <= gamma) throw InvalidArgumentException(HERE) << "mu must be greater than gamma, here mu=" << mu << " and gamma=" << gamma;
108
109 Point muSigmaParametersValues(inP);
110 muSigmaParametersValues[1] *= mu;
111 const LogNormalMuSigma muSigmaParameters(mu, sigmaOverMu * mu, gamma);
112
113 return muSigmaParameters(muSigmaParametersValues);
114 }
115
116
inverse(const Point & inP) const117 Point LogNormalMuSigmaOverMu::inverse(const Point & inP) const
118 {
119 const LogNormalMuSigma muSigmaParameters;
120 Point muSigmaOverMuParameters(muSigmaParameters.inverse(inP));
121 const Scalar mu = muSigmaOverMuParameters[0];
122 if (mu == 0.0) throw InvalidArgumentException(HERE) << "Error: mu cannot be null in the parameter set (mu, sigmaOverMu)";
123 muSigmaOverMuParameters[1] /= mu;
124 return muSigmaOverMuParameters;
125 }
126
127 /* Parameters value and description accessor */
setValues(const Point & inP)128 void LogNormalMuSigmaOverMu::setValues(const Point & inP)
129 {
130 if (inP.getDimension() != 3) throw InvalidArgumentException(HERE) << "the given point must have dimension=3, here dimension=" << inP.getDimension();
131 mu_ = inP[0];
132 sigmaOverMu_ = inP[1];
133 gamma_ = inP[2];
134 }
135
getValues() const136 Point LogNormalMuSigmaOverMu::getValues() const
137 {
138 Point point(3);
139 point[0] = mu_;
140 point[1] = sigmaOverMu_;
141 point[2] = gamma_;
142 return point;
143 }
144
getDescription() const145 Description LogNormalMuSigmaOverMu::getDescription() const
146 {
147 Description description(3);
148 description[0] = "mu";
149 description[1] = "sigmaOverMu";
150 description[2] = "gamma";
151 return description;
152 }
153
154 /* String converter */
__repr__() const155 String LogNormalMuSigmaOverMu::__repr__() const
156 {
157 OSS oss(true);
158 oss << "class=" << GetClassName()
159 << " name=" << getName()
160 << " mu=" << mu_
161 << " sigmaOverMu=" << sigmaOverMu_
162 << " gamma=" << gamma_;
163 return oss;
164 }
165
166
__str__(const String &) const167 String LogNormalMuSigmaOverMu::__str__(const String & ) const
168 {
169 OSS oss(false);
170 oss << getClassName() << "(mu = " << mu_ << ", sigmaOverMu = " << sigmaOverMu_ << ", gamma = " << gamma_ << ")";
171 return oss;
172 }
173
save(Advocate & adv) const174 void LogNormalMuSigmaOverMu::save(Advocate & adv) const
175 {
176 DistributionParametersImplementation::save(adv);
177 adv.saveAttribute( "mu_", mu_ );
178 adv.saveAttribute( "sigmaOverMu_", sigmaOverMu_ );
179 adv.saveAttribute( "gamma_", gamma_ );
180 }
181
load(Advocate & adv)182 void LogNormalMuSigmaOverMu::load(Advocate & adv)
183 {
184 DistributionParametersImplementation::load(adv);
185 adv.loadAttribute( "mu_", mu_ );
186 adv.loadAttribute( "sigmaOverMu_", sigmaOverMu_ );
187 adv.loadAttribute( "gamma_", gamma_ );
188 }
189
190 END_NAMESPACE_OPENTURNS
191