Skip to content

Commit d813768

Browse files
IthanilIthanil
authored
Merge pull request #35 from DCM-UPB/fix_euclidmaps
Fix euclidmaps
2 parents 848df1f + bb4989b commit d813768

File tree

4 files changed

+143
-105
lines changed

4 files changed

+143
-105
lines changed

‎src/feature_maps/EuclideanDistanceMap.cpp‎

Lines changed: 55 additions & 49 deletions
Original file line numberDiff line numberDiff line change
@@ -15,7 +15,7 @@ double calcDistHelper(const std::vector<double> &srcv, const std::vector<double>
1515
for (size_t i=0; i<ndim; ++i) {
1616
dist += pow(srcv[i] - fixedPoint[i], 2);
1717
}
18-
return sqrt(dist);
18+
return dist;
1919
}
2020

2121
double EuclideanDistanceMap::_calcDist()
@@ -24,7 +24,7 @@ double EuclideanDistanceMap::_calcDist()
2424
for (size_t i=0; i<_ndim; ++i) {
2525
dist += pow(_sources[i]->getValue() - _fixedPoint[i], 2);
2626
}
27-
return sqrt(dist);
27+
return dist;
2828
}
2929

3030

@@ -87,14 +87,13 @@ double EuclideanDistanceMap::getFeedSigma()
8787
for (size_t i=0; i<_sources.size(); ++i) {
8888
srcv.push_back(_sources[i]->getOutputMu());
8989
}
90-
double dist = calcDistHelper(srcv, _fixedPoint, _ndim);
9190

9291
double sigma = 0.;
9392
for (size_t i=0; i<_ndim; ++i) {
94-
sigma += srcv[i]*srcv[i] * pow(_sources[i]->getOutputSigma(), 2); // (d²/dx² sigmaX)²
93+
sigma += 4.0 * srcv[i]*srcv[i] * pow(_sources[i]->getOutputSigma(), 2); // (d²/dx² sigmaX)²
9594
}
9695

97-
return sqrt(sigma) / dist;
96+
return sqrt(sigma);
9897
}
9998

10099

@@ -109,80 +108,87 @@ double EuclideanDistanceMap::getFeed()
109108

110109
double EuclideanDistanceMap::getFirstDerivativeFeed(const int &i1d)
111110
{
112-
const double dist = _calcDist();
113-
114-
if (dist > 0) {
115-
double d1 = 0.;
116-
for (size_t i=0; i<_ndim; ++i) {
117-
const double v = _sources[i]->getValue();
118-
const double d1v = _sources[i]->getFirstDerivativeValue(i1d);
119-
d1 += (v - _fixedPoint[i]) * d1v;
120-
}
121-
122-
return d1 / dist;
111+
double d1 = 0.;
112+
for (size_t i=0; i<_ndim; ++i) {
113+
const double v = _sources[i]->getValue();
114+
const double d1v = _sources[i]->getFirstDerivativeValue(i1d);
115+
d1 += (v - _fixedPoint[i]) * d1v;
123116
}
124-
else return 0.;
117+
118+
return 2.0 * d1;
125119
}
126120

127121

128122
double EuclideanDistanceMap::getSecondDerivativeFeed(const int &i2d)
129123
{
130-
const double dist = _calcDist();
131-
132-
if (dist > 0) {
133-
const double dist2 = dist * dist;
134-
135-
double d21 = 0.;
136-
double d22 = 0.;
137-
for (size_t i=0; i<_ndim; ++i) {
138-
const double v = _sources[i]->getValue();
139-
const double d1v = _sources[i]->getFirstDerivativeValue(i2d);
140-
const double d2v = _sources[i]->getSecondDerivativeValue(i2d);
141-
142-
const double dv = v - _fixedPoint[i];
143-
d21 += dv * d2v + d1v*d1v;
144-
d22 += dv * d1v;
145-
}
146-
147-
return (d21 - d22*d22 / dist2) / dist;
124+
double d2 = 0.;
125+
for (size_t i=0; i<_ndim; ++i) {
126+
const double v = _sources[i]->getValue();
127+
const double d1v = _sources[i]->getFirstDerivativeValue(i2d);
128+
const double d2v = _sources[i]->getSecondDerivativeValue(i2d);
129+
d2 += d1v*d1v + (v - _fixedPoint[i]) * d2v;
148130
}
149-
else return 0.;
131+
132+
return 2.0 * d2;
150133
}
151134

135+
152136
double EuclideanDistanceMap::getVariationalFirstDerivativeFeed(const int &iv1d)
137+
153138
{
154139
if (iv1d < _vp_id_shift) {
155-
const double dist = _calcDist();
156-
157-
if (dist > 0) {
158-
double vd1 = 0.;
159-
for (size_t i=0; i<_ndim; ++i) {
160-
const double v = _sources[i]->getValue();
161-
const double vdv = _sources[i]->getVariationalFirstDerivativeValue(iv1d);
162-
163-
vd1 += (v - _fixedPoint[i]) * vdv;
164-
}
140+
double vd1 = 0.;
141+
for (size_t i=0; i<_ndim; ++i) {
142+
const double v = _sources[i]->getValue();
143+
const double vdv = _sources[i]->getVariationalFirstDerivativeValue(iv1d);
165144

166-
return vd1 / dist;
145+
vd1 += (v - _fixedPoint[i]) * vdv;
167146
}
147+
148+
return 2.0 * vd1;
168149
}
150+
169151
return 0.;
170152
}
171153

172154

173155
double EuclideanDistanceMap::getCrossFirstDerivativeFeed(const int &i1d, const int &iv1d)
174156
{
175157
if (iv1d < _vp_id_shift) {
176-
throw std::runtime_error("[EuclideanDistanceMap::getCrossFirstDerivativeFeed] Variational cross derivatives not supported yet.");
158+
double cd1 = 0.;
159+
for (size_t i=0; i<_ndim; ++i) {
160+
const double v = _sources[i]->getValue();
161+
const double d1v = _sources[i]->getFirstDerivativeValue(i1d);
162+
const double vdv = _sources[i]->getVariationalFirstDerivativeValue(iv1d);
163+
const double cd1v = _sources[i]->getCrossFirstDerivativeValue(i1d, iv1d);
164+
165+
cd1 += d1v * vdv + (v - _fixedPoint[i]) * cd1v;
166+
}
167+
168+
return 2.0 * cd1;
177169
}
170+
178171
else return 0.;
179172
}
180173

181174

182175
double EuclideanDistanceMap::getCrossSecondDerivativeFeed(const int &i2d, const int &iv2d)
183176
{
184177
if (iv2d < _vp_id_shift) {
185-
throw std::runtime_error("[EuclideanDistanceMap::getCrossFirstDerivativeFeed] Variational cross derivatives not supported yet.");
178+
double cd2 = 0.;
179+
for (size_t i=0; i<_ndim; ++i) {
180+
const double v = _sources[i]->getValue();
181+
const double d1v = _sources[i]->getFirstDerivativeValue(i2d);
182+
const double d2v = _sources[i]->getSecondDerivativeValue(i2d);
183+
const double vdv = _sources[i]->getVariationalFirstDerivativeValue(iv2d);
184+
const double cd1v = _sources[i]->getCrossFirstDerivativeValue(i2d, iv2d);
185+
const double cd2v = _sources[i]->getCrossFirstDerivativeValue(i2d, iv2d);
186+
187+
cd2 += d2v * vdv + d1v * cd1v + d1v * cd1v + (v - _fixedPoint[i]) * cd2v;
188+
}
189+
190+
return 2.0 * cd2;
186191
}
192+
187193
else return 0.;
188194
}

‎src/feature_maps/EuclideanDistanceMap.hpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -7,6 +7,7 @@
77
#include <vector>
88
#include <cstddef>
99

10+
// takes coordinates on input side and calculates squared euclidean distance against a fixed point
1011
class EuclideanDistanceMap: public MultiDimStaticMap
1112
{
1213
protected:

‎src/feature_maps/EuclideanPairDistanceMap.cpp‎

Lines changed: 86 additions & 56 deletions
Original file line numberDiff line numberDiff line change
@@ -15,7 +15,7 @@ double calcDistHelper(std::vector<double> &srcv, const size_t &ndim)
1515
for (size_t i=0; i<ndim; ++i) {
1616
dist += pow(srcv[i] - srcv[i+ndim], 2);
1717
}
18-
return sqrt(dist);
18+
return dist;
1919
}
2020

2121
double EuclideanPairDistanceMap::_calcDist()
@@ -24,7 +24,7 @@ double EuclideanPairDistanceMap::_calcDist()
2424
for (size_t i=0; i<_ndim; ++i) {
2525
dist += pow(_sources[i]->getValue() - _sources[i+_ndim]->getValue(), 2);
2626
}
27-
return sqrt(dist);
27+
return dist;
2828
}
2929

3030

@@ -55,14 +55,13 @@ double EuclideanPairDistanceMap::getFeedSigma()
5555
for (size_t i=0; i<_sources.size(); ++i) {
5656
srcv.push_back(_sources[i]->getOutputMu());
5757
}
58-
double dist = calcDistHelper(srcv, _ndim);
5958

6059
double sigma = 0.;
6160
for (size_t i=0; i<_ndim; ++i) {
62-
sigma += pow(srcv[i] - srcv[i+_ndim], 2) * ( pow(_sources[i]->getOutputSigma(), 2) + pow(_sources[i+_ndim]->getOutputSigma(), 2) ); // (d²/dx² sigmaX)²
61+
sigma += 2.0 * pow(srcv[i] - srcv[i+_ndim], 2) * ( pow(_sources[i]->getOutputSigma(), 2) + pow(_sources[i+_ndim]->getOutputSigma(), 2) ); // (d²/dx² sigmaX)²
6362
}
6463

65-
return sqrt(sigma) / dist;
64+
return sqrt(sigma);
6665
}
6766

6867

@@ -77,88 +76,119 @@ double EuclideanPairDistanceMap::getFeed()
7776

7877
double EuclideanPairDistanceMap::getFirstDerivativeFeed(const int &i1d)
7978
{
80-
const double dist = _calcDist();
81-
82-
if (dist > 0) {
83-
double d1 = 0.;
84-
for (size_t i=0; i<_ndim; ++i) {
85-
const double v1 = _sources[i]->getValue();
86-
const double v2 = _sources[i+_ndim]->getValue();
87-
const double d1v1 = _sources[i]->getFirstDerivativeValue(i1d);
88-
const double d1v2 = _sources[i+_ndim]->getFirstDerivativeValue(i1d);
89-
d1 += (v1 - v2) * (d1v1 - d1v2);
90-
}
91-
92-
return d1 / dist;
79+
double d1 = 0.;
80+
for (size_t i=0; i<_ndim; ++i) {
81+
const double v1 = _sources[i]->getValue();
82+
const double v2 = _sources[i+_ndim]->getValue();
83+
const double d1v1 = _sources[i]->getFirstDerivativeValue(i1d);
84+
const double d1v2 = _sources[i+_ndim]->getFirstDerivativeValue(i1d);
85+
d1 += (v1 - v2) * (d1v1 - d1v2);
9386
}
94-
else return 0.;
87+
88+
return 2.0 * d1;
9589
}
9690

9791

9892
double EuclideanPairDistanceMap::getSecondDerivativeFeed(const int &i2d)
9993
{
100-
const double dist = _calcDist();
94+
double d2 = 0.;
95+
for (size_t i=0; i<_ndim; ++i) {
96+
const double v1 = _sources[i]->getValue();
97+
const double v2 = _sources[i+_ndim]->getValue();
98+
const double d1v1 = _sources[i]->getFirstDerivativeValue(i2d);
99+
const double d1v2 = _sources[i+_ndim]->getFirstDerivativeValue(i2d);
100+
const double d2v1 = _sources[i]->getSecondDerivativeValue(i2d);
101+
const double d2v2 = _sources[i+_ndim]->getSecondDerivativeValue(i2d);
102+
103+
const double diff = v1 - v2;
104+
const double d1diff = d1v1 - d1v2;
105+
const double d2diff = d2v1 - d2v2;
106+
107+
d2 += d1diff * d1diff + diff * d2diff;
108+
}
101109

102-
if (dist > 0) {
103-
const double dist2 = dist * dist;
110+
return 2.0 * d2;
111+
}
104112

105-
double d21 = 0.;
106-
double d22 = 0.;
113+
double EuclideanPairDistanceMap::getVariationalFirstDerivativeFeed(const int &iv1d)
114+
{
115+
if (iv1d < _vp_id_shift) {
116+
double vd1 = 0.;
107117
for (size_t i=0; i<_ndim; ++i) {
108118
const double v1 = _sources[i]->getValue();
109119
const double v2 = _sources[i+_ndim]->getValue();
110-
const double d1v1 = _sources[i]->getFirstDerivativeValue(i2d);
111-
const double d1v2 = _sources[i+_ndim]->getFirstDerivativeValue(i2d);
112-
const double d2v1 = _sources[i]->getSecondDerivativeValue(i2d);
113-
const double d2v2 = _sources[i+_ndim]->getSecondDerivativeValue(i2d);
120+
const double vdv1 = _sources[i]->getVariationalFirstDerivativeValue(iv1d);
121+
const double vdv2 = _sources[i+_ndim]->getVariationalFirstDerivativeValue(iv1d);
114122

115-
const double dv = (v1 -v2);
116-
const double d1d = d1v1 - d1v2;
117-
const double d2d = d2v1 - d2v2;
118-
d21 += dv * d2d + d1d*d1d;
119-
d22 += dv * d1d;
123+
vd1 += (v1 - v2) * (vdv1 - vdv2);
120124
}
121125

122-
return (d21 - d22*d22 / dist2) / dist;
126+
return 2.0 * vd1;
123127
}
124-
else return 0.;
125-
}
126128

127-
double EuclideanPairDistanceMap::getVariationalFirstDerivativeFeed(const int &iv1d)
128-
{
129-
if (iv1d < _vp_id_shift) {
130-
const double dist = _calcDist();
131-
132-
if (dist > 0) {
133-
double vd1 = 0.;
134-
for (size_t i=0; i<_ndim; ++i) {
135-
const double v1 = _sources[i]->getValue();
136-
const double v2 = _sources[i+_ndim]->getValue();
137-
const double vdv1 = _sources[i]->getVariationalFirstDerivativeValue(iv1d);
138-
const double vdv2 = _sources[i+_ndim]->getVariationalFirstDerivativeValue(iv1d);
139-
vd1 += (v1 - v2) * (vdv1-vdv2);
140-
}
141-
142-
return vd1 / dist;
143-
}
144-
}
145129
return 0.;
146130
}
147131

148132

149133
double EuclideanPairDistanceMap::getCrossFirstDerivativeFeed(const int &i1d, const int &iv1d)
150134
{
151135
if (iv1d < _vp_id_shift) {
152-
throw std::runtime_error("[EuclideanPairDistanceMap::getCrossFirstDerivativeFeed] Variational cross derivatives not supported yet.");
136+
double cd1 = 0.;
137+
for (size_t i=0; i<_ndim; ++i) {
138+
const double v1 = _sources[i]->getValue();
139+
const double v2 = _sources[i+_ndim]->getValue();
140+
const double d1v1 = _sources[i]->getFirstDerivativeValue(i1d);
141+
const double d1v2 = _sources[i+_ndim]->getFirstDerivativeValue(i1d);
142+
const double vdv1 = _sources[i]->getVariationalFirstDerivativeValue(iv1d);
143+
const double vdv2 = _sources[i+_ndim]->getVariationalFirstDerivativeValue(iv1d);
144+
const double cd1v1 = _sources[i]->getCrossFirstDerivativeValue(i1d, iv1d);
145+
const double cd1v2 = _sources[i+_ndim]->getCrossFirstDerivativeValue(i1d, iv1d);
146+
147+
const double diff = v1 - v2;
148+
const double d1diff = d1v1 - d1v2;
149+
const double vddiff = vdv1 - vdv2;
150+
const double cd1diff = cd1v1 - cd1v2;
151+
152+
cd1 += d1diff * vddiff + diff * cd1diff;
153+
}
154+
155+
return 2.0 * cd1;
153156
}
157+
154158
else return 0.;
155159
}
156160

157161

158162
double EuclideanPairDistanceMap::getCrossSecondDerivativeFeed(const int &i2d, const int &iv2d)
159163
{
160164
if (iv2d < _vp_id_shift) {
161-
throw std::runtime_error("[EuclideanPairDistanceMap::getCrossFirstDerivativeFeed] Variational cross derivatives not supported yet.");
165+
double cd2 = 0.;
166+
for (size_t i=0; i<_ndim; ++i) {
167+
const double v1 = _sources[i]->getValue();
168+
const double v2 = _sources[i+_ndim]->getValue();
169+
const double d1v1 = _sources[i]->getFirstDerivativeValue(i2d);
170+
const double d1v2 = _sources[i+_ndim]->getFirstDerivativeValue(i2d);
171+
const double d2v1 = _sources[i]->getSecondDerivativeValue(i2d);
172+
const double d2v2 = _sources[i+_ndim]->getSecondDerivativeValue(i2d);
173+
const double vdv1 = _sources[i]->getVariationalFirstDerivativeValue(iv2d);
174+
const double vdv2 = _sources[i+_ndim]->getVariationalFirstDerivativeValue(iv2d);
175+
const double cd1v1 = _sources[i]->getCrossFirstDerivativeValue(i2d, iv2d);
176+
const double cd1v2 = _sources[i+_ndim]->getCrossFirstDerivativeValue(i2d, iv2d);
177+
const double cd2v1 = _sources[i]->getCrossSecondDerivativeValue(i2d, iv2d);
178+
const double cd2v2 = _sources[i+_ndim]->getCrossSecondDerivativeValue(i2d, iv2d);
179+
180+
const double diff = v1 - v2;
181+
const double d1diff = d1v1 - d1v2;
182+
const double d2diff = d2v1 - d2v2;
183+
const double vddiff = vdv1 - vdv2;
184+
const double cd1diff = cd1v1 - cd1v2;
185+
const double cd2diff = cd2v1 - cd2v2;
186+
187+
cd2 += d2diff * vddiff + 2.0 * d1diff * cd1diff + diff * cd2diff;
188+
}
189+
190+
return 2.0 * cd2;
162191
}
192+
163193
else return 0.;
164194
}

‎src/feature_maps/EuclideanPairDistanceMap.hpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -4,6 +4,7 @@
44
#include "MultiDimStaticMap.hpp"
55
#include "NetworkLayer.hpp"
66

7+
// takes coordinates on input side and calculates squared euclidean distance of a pair
78
class EuclideanPairDistanceMap: public MultiDimStaticMap
89
{
910
protected:

0 commit comments

Comments
 (0)