Commit 9f384106 authored by incardon's avatar incardon

Adding PS-CMA-ES

parent 8c7a988c
master
This diff is collapsed.
/*
* f15_cec_fun.hpp
*
* Created on: Jan 14, 2018
* Author: i-bird
*/
#ifndef EXAMPLE_NUMERICS_PS_CMA_ES_F15_CEC_FUN_HPP_
#define EXAMPLE_NUMERICS_PS_CMA_ES_F15_CEC_FUN_HPP_
#include "f15_cec_const.hpp"
#include <limits>
#include <math.h>
template<unsigned int dim>
void Job15(int funcnr,Eigen::VectorXd & vars,double & res)
{
// local used vars
double sum,sum1,sum2,prod,e1,e2;
int i,j,k;
// weierstrass vars
int Kmax = 20;
const double a_c = 0.5;
const double b_c = 3.0;
if (funcnr < 2)
{
// rastrigin
sum = 10.0 * dim;
for (size_t i = 0 ; i < dim ; i++)
{
sum += vars(i)*vars(i);
sum -= 10.0*cos(2*M_PI*vars[i]);
}
res = sum;
}
else if (funcnr < 4)
{
// weierstrass
sum1 = 0.0;
sum2 = 0.0;
double a_k = 1.0;
double b_k = 1.0;
for (size_t i = 0 ; i < dim ; i++)
{
a_k = 1.0;
b_k = 1.0;
for (size_t j = 0 ; j <= Kmax ; j++, a_k *= a_c,b_k *= b_c)
{
sum1 = sum1 + a_k * cos((M_PI)*2.0 * b_k * (vars(i)+0.5));
}
}
a_k = 1.0;
b_k = 1.0;
for (size_t j = 0 ; j <= Kmax ; j++, a_k *= a_c, b_k *= b_c)
{
sum2 = sum2 + a_k * cos((M_PI)*2.0 * b_k * (0.5));
}
res = sum1 - sum2*dim;
}
else if (funcnr < 6)
{
// griewank
prod = 1;
sum = 0.0;
for (size_t i = 1 ; i <= dim ; i++)
{
sum= sum + (vars(i-1)*vars(i-1))/4000.0;
prod=prod * cos(vars(i-1)/(sqrt(double(i))));
}
res = sum-prod+1;
}
else if (funcnr < 8)
{
// ackley
e1 = 0.0;
e2 = 0.0;
for (size_t i = 0 ; i < dim ; i++)
{
e1 = e1 + vars(i)*vars(i);
e2 = e2 + cos(2.0*M_PI*vars(i));
}
res = exp(1.0) + 20.0 - 20*exp(-0.2*sqrt(e1/dim));
res = res - exp(e2/dim);
}
else if (funcnr <= 10)
{
// sphere
sum = vars.transpose() * vars;
res = sum;
}
}
template<unsigned int dim>
double hybrid_composition(Eigen::VectorXd & vars)
{
double ZBQLNOR;
//local used vars
double wMax,sumSqr,wSum,w1mMaxPow;
int i,j,k;
double sumF,t_res;
Eigen::VectorXd job_z[10];
for (size_t i = 0 ; i < 10 ; i++)
{job_z[i].resize(dim);}
double job_w[10];
double res = 0.0;
for (size_t i = 0 ; i < dim ; i++)
{
if (vars[i] < -5.0 || vars[i] > 5.0)
{return std::numeric_limits<double>::infinity();}
}
// get the raw weights
wMax = - std::numeric_limits<double>::max();
for (size_t i = 0; i < 10 ; i++)
{
sumSqr = 0.0;
//Shift the Input
job_z[i] = vars - f15_o[i];
sumSqr += (job_z[i].transpose() * job_z[i]);
job_w[i] = exp(-1.0 * sumSqr / (2.0 * dim));
if (wMax < job_w[i])
{wMax = job_w[i];}
}
// Modify the weights
wSum = 0.0;
w1mMaxPow = 1.0 - wMax*wMax*wMax*wMax*wMax*wMax*wMax*wMax*wMax*wMax;
for (size_t i = 0; i < 10 ; i++)
{
if (job_w[i] != wMax)
{job_w[i] = job_w[i]* w1mMaxPow;};
wSum = wSum + job_w[i];
}
// Normalize the weights
for (size_t i = 0; i < 10 ; i++)
{job_w[i] /= wSum;}
sumF = 0.0;
for (size_t i = 0; i < 10 ; i++)
{
job_z[i] = job_z[i] / job_lambda[i];
//calling the basic functions
Job15<dim>(i,job_z[i],t_res);
sumF = sumF + job_w[i] * (2000.0*t_res/f15_max[i] + bias[i]);
}
res = sumF + 120;
return res;
}
template<unsigned int dim>
void prepare_f15()
{
// load f15_o
for (size_t j = 0 ; j < 10 ; j++)
{
Eigen::VectorXd fmp(dim);
f15_o[j].resize(dim);
for (size_t i = 0 ; i < dim ; i++)
{
f15_o[j](i) = f15_const[j][i];
fmp(i) = 5.0 / job_lambda[j];
}
double result;
Job15<dim>(j,fmp,result);
f15_max[j] = fabs(result);
}
}
#endif /* EXAMPLE_NUMERICS_PS_CMA_ES_F15_CEC_FUN_HPP_ */
This diff is collapsed.
......@@ -877,9 +877,9 @@ template<typename grid> void calc_rhs(grid & g_vort, grid & g_vel, grid & g_dwp)
// calculate several pre-factors for the stencil finite
// difference
float fac1 = 2.0f*nu/(g_vort.spacing(0)*g_vort.spacing(0));
float fac2 = 2.0f*nu/(g_vort.spacing(1)*g_vort.spacing(1));
float fac3 = 2.0f*nu/(g_vort.spacing(2)*g_vort.spacing(2));
float fac1 = 1.0f*nu/(g_vort.spacing(0)*g_vort.spacing(0));
float fac2 = 1.0f*nu/(g_vort.spacing(1)*g_vort.spacing(1));
float fac3 = 1.0f*nu/(g_vort.spacing(2)*g_vort.spacing(2));
float fac4 = 0.5f/(g_vort.spacing(0));
float fac5 = 0.5f/(g_vort.spacing(1));
......
......@@ -121,7 +121,7 @@ class domain_nn_calculator_cart
// +2 is padding
for (size_t j = 0 ; j < dim ; j++)
sz[j] = proc_box.getHigh(j) - proc_box.getLow(j) + 2 + 1;
{sz[j] = proc_box.getHigh(j) - proc_box.getLow(j) + 2 + 1;}
gs.setDimensions(sz);
......@@ -130,7 +130,7 @@ class domain_nn_calculator_cart
g.setMemory();
for (size_t i = 0 ; i < dim ; i++)
one.set_d(i,1);
{one.set_d(i,1);}
// Calculate the csr neighborhood
openfpm::vector<std::pair<grid_key_dx<dim>,grid_key_dx<dim>>> csr;
......@@ -191,7 +191,7 @@ class domain_nn_calculator_cart
sub_keys.last().NN_subsub.resize(g.template get<0>(key).size());
for (size_t i = 0 ; i < g.template get<0>(key).size() ; i++)
sub_keys.last().NN_subsub.get(i) = g.template get<0>(key).get(i) - one;
{sub_keys.last().NN_subsub.get(i) = g.template get<0>(key).get(i) - one;}
}
++it;
......
......@@ -83,6 +83,12 @@ struct Box_sub
//! see ie_ghost follow sector explanation
comb<dim> cmb;
//! Constructor reset cmb
Box_sub()
{
cmb.zero();
}
};
//! Particular case for local internal ghost boxes
......
......@@ -602,18 +602,6 @@ class grid_dist_id : public grid_dist_id_comm<dim,St,T,Decomposition,Memory,devi
protected:
/*! \brief Get the point where it start the origin of the grid of the sub-domain i
*
* \param i sub-domain
*
* \return the point
*
*/
Point<dim,St> getOffset(size_t i)
{
return pmul(Point<dim,St>(gdb_ext.get(i).origin), cd_sm.getCellBox().getP2());
}
/*! \brief Given a local sub-domain i with a local grid Domain + ghost return the part of the local grid that is domain
*
* \param i sub-domain
......@@ -687,6 +675,18 @@ public:
return domain;
}
/*! \brief Get the point where it start the origin of the grid of the sub-domain i
*
* \param i sub-domain
*
* \return the point
*
*/
Point<dim,St> getOffset(size_t i)
{
return pmul(Point<dim,St>(gdb_ext.get(i).origin), cd_sm.getCellBox().getP2()) + getDomain().getP1();
}
/*! \brief Get the spacing of the grid in direction i
*
* \param i dimension
......
......@@ -568,11 +568,6 @@ public:
while (it.isNext())
{
//auto key = it.get();
//if (g.template get<0>(key) != 1)
//std::cout << "WRONG???????" << std::endl;
++it;
count++;
}
......@@ -581,13 +576,11 @@ public:
Point<dim,St> p;
for (size_t n = 0; n < dim; n++)
p.get(n) = g.getGrid().getBox().getHigh(n);
//std::cout << "G after resize: (" << g.getGrid().getBox().getLow(0) << "; " << g.getGrid().getBox().getLow(1) << "); (" << g.getGrid().getBox().getHigh(0) << "; " << g.getGrid().getBox().getHigh(1) << ")" << std::endl;
{p.get(n) = g.getGrid().getBox().getHigh(n);}
Point<dim,St> point;
for (size_t n = 0; n < dim; n++)
point.get(n) = (b.getHigh(n) + b.getLow(n))/2;
{point.get(n) = (b.getHigh(n) + b.getLow(n))/2;}
for (size_t j = 0; j < gdb_ext.size(); j++)
{
......@@ -612,7 +605,6 @@ public:
std::string str = key.to_string();
grid_key_dx<dim> key2 = key - start;
//std::cout << "Key: " << str << std::endl;
loc_grid.get(j).get_o(key) = g.get_o(key2);
count2++;
......@@ -622,7 +614,6 @@ public:
}
}
}
//std::cout << "Count after: " << count2 << std::endl;
}
/*! \brief Label intersection grids for mappings
......
......@@ -1871,6 +1871,27 @@ BOOST_AUTO_TEST_CASE ( grid_ghost_correction )
Test_ghost_correction(domain,k,4);
}
BOOST_AUTO_TEST_CASE ( grid_basic_functions )
{
auto & v_cl = create_vcluster();
if (v_cl.getProcessingUnits() != 1)
{return;}
size_t sz[2] = {(size_t)8,(size_t)8};
periodicity<2> bc = {PERIODIC,PERIODIC};
Ghost<2,long int> g(1);
Box<2,double> domain({-1.0,-1.0},{1.0,1.0});
grid_dist_id<2, double, aggregate<double>> grid(sz,domain,g,bc);
std::cout << "Offset: " << grid.getOffset(0).toString() << std::endl;
BOOST_REQUIRE_EQUAL(grid.getOffset(0)[0],-1.25);
BOOST_REQUIRE_EQUAL(grid.getOffset(0)[1],-1.25);
}
BOOST_AUTO_TEST_CASE ( grid_overflow_round_off_error )
{
size_t numGridPoint = 100;
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment