mskefile，データ例，プログラム
-------------------------i_data-------------------------
// 関数 a
関数 1 変数の数 2 最大試行回数 1000 初期値設定係数 1.0
許容誤差 1.0e-20 縮小比率 0.7 拡大比率 1.5
初期値 0.0 0.0
// 関数 b
関数 2 変数の数 2 最大試行回数 1000 初期値設定係数 1.0
許容誤差 1.0e-20 縮小比率 0.7 拡大比率 1.5
初期値 0.0 0.0
// 関数 c
関数 3 変数の数 2 最大試行回数 1000 初期値設定係数 1.0
許容誤差 1.0e-20 縮小比率 0.7 拡大比率 1.5
初期値 0.0 0.0
-------------------------program-------------------------
/**************************************/
/* シンプレックス法による最小値の計算 */
/* coded by Y.Suganuma */
/**************************************/
#include <stdio.h>
double snx1(double *);
double snx2(double *);
double snx3(double *);
int simplex(int, int, double, double, double, double, double *, double *, double (*)(double *));
int main()
{
double eps, r1, r2, *x, y, k;
int fun, i1, max, n, sw = 0;
// データの入力
scanf("%*s %d %*s %d %*s %d %*s %lf", &fun, &n, &max, &k);
scanf("%*s %lf %*s %lf %*s %lf", &eps, &r1, &r2);
x = new double [n];
scanf("%*s");
for (i1 = 0; i1 < n; i1++)
scanf("%lf", &x[i1]);
// 実行
switch (fun) {
case 1:
sw = simplex(max, n, k, eps, r1, r2, &y, x, snx1);
break;
case 2:
sw = simplex(max, n, k, eps, r1, r2, &y, x, snx2);
break;
case 3:
sw = simplex(max, n, k, eps, r1, r2, &y, x, snx3);
break;
}
// 結果の出力
if (sw < 0)
printf(" 収束しませんでした！");
else {
printf(" 結果＝");
for (i1 = 0; i1 < n; i1++)
printf("%f ", x[i1]);
printf(" 最小値＝%f 回数＝%d\n", y, sw);
}
delete [] x;
return 0;
}
/******************************************/
/* シンプレックス法 */
/* max : 最大繰り返し回数 */
/* n : 次元 */
/* k : 初期値設定係数 */
/* r1 : 縮小比率 */
/* r2 : 拡大比率 */
/* y : 最小値 */
/* x : x(初期値と答え) */
/* snx : 関数値を計算する関数名 */
/* return : >=0 : 正常終了(収束回数) */
/* =-1 : 収束せず */
/******************************************/
#include <math.h>
int simplex(int max, int n, double k, double eps, double r1, double r2, double *y, double *x, double (*snx)(double *))
{
double **xx, *yy, *xg, *xr, *xn, yr, yn, e;
int i1, i2, fh = -1, fg = -1, fl = -1, count = 0, sw = -1;
// 初期値の設定
yy = new double [n+1];
xg = new double [n];
xr = new double [n];
xn = new double [n];
xx = new double * [n+1];
for (i1 = 0; i1 < n+1; i1++)
xx[i1] = new double [n];
for (i1 = 0; i1 < n+1; i1++) {
for (i2 = 0; i2 < n; i2++)
xx[i1][i2] = x[i2];
if (i1 > 0)
xx[i1][i1-1] += k;
yy[i1] = snx(xx[i1]);
}
// 最大値，最小値の計算
for (i1 = 0; i1 < n+1; i1++) {
if (fh < 0 || yy[i1] > yy[fh])
fh = i1;
if (fl < 0 || yy[i1] < yy[fl])
fl = i1;
}
for (i1 = 0; i1 < n+1; i1++) {
if (i1 != fh && (fg < 0 || yy[i1] > yy[fg]))
fg = i1;
}
// 最小値の計算
while (count < max && sw < 0) {
count++;
// 重心の計算
for (i1 = 0; i1 < n; i1++)
xg[i1] = 0.0;
for (i1 = 0; i1 < n+1; i1++) {
if (i1 != fh) {
for (i2 = 0; i2 < n; i2++)
xg[i2] += xx[i1][i2];
}
}
for (i1 = 0; i1 < n; i1++)
xg[i1] /= n;
// 最大点の置き換え
for (i1 = 0; i1 < n; i1++)
xr[i1] = 2.0 * xg[i1] - xx[fh][i1];
yr = snx(xr);
if (yr >= yy[fh]) { // 縮小
for (i1 = 0; i1 < n; i1++)
xr[i1] = (1.0 - r1) * xx[fh][i1] + r1 * xr[i1];
yr = snx(xr);
}
else if (yr < (yy[fl]+(r2-1.0)*yy[fh])/r2) { // 拡大
for (i1 = 0; i1 < n; i1++)
xn[i1] = r2 * xr[i1] - (r2 - 1.0) * xx[fh][i1];
yn = snx(xn);
if (yn <= yr) {
for (i1 = 0; i1 < n; i1++)
xr[i1] = xn[i1];
yr = yn;
}
}
for (i1 = 0; i1 < n; i1++)
xx[fh][i1] = xr[i1];
yy[fh] = yr;
// シンプレックス全体の縮小
if (yy[fh] >= yy[fg]) {
for (i1 = 0; i1 < n+1; i1++) {
if (i1 != fl) {
for (i2 = 0; i2 < n; i2++)
xx[i1][i2] = 0.5 * (xx[i1][i2] + xx[fl][i2]);
yy[i1] = snx(xx[i1]);
}
}
}
// 最大値，最小値の計算
fh = -1;
fg = -1;
fl = -1;
for (i1 = 0; i1 < n+1; i1++) {
if (fh < 0 || yy[i1] > yy[fh])
fh = i1;
if (fl < 0 || yy[i1] < yy[fl])
fl = i1;
}
for (i1 = 0; i1 < n+1; i1++) {
if (i1 != fh && (fg < 0 || yy[i1] > yy[fg]))
fg = i1;
}
// 収束判定
e = 0.0;
for (i1 = 0; i1 < n+1; i1++) {
if (i1 != fl) {
yr = yy[i1] - yy[fl];
e += yr * yr;
}
}
if (e < eps)
sw = 0;
}
if (sw == 0) {
sw = count;
*y = yy[fl];
for (i1 = 0; i1 < n; i1++)
x[i1] = xx[fl][i1];
}
delete [] yy;
delete [] xg;
delete [] xr;
delete [] xn;
for (i1 = 0; i1 < n+1; i1++)
delete [] xx[i1];
delete [] xx;
return sw;
}
/*******************/
/* 関数値の計算 */
/* y = f(x) */
/* return : y */
/*******************/
// 関数１
double snx1(double *x)
{
double x1, y1, y;
x1 = x[0] - 1.0;
y1 = x[1] - 2.0;
y = x1 * x1 + y1 * y1;
return y;
}
// 関数２
double snx2(double *x)
{
double x1, y1, y;
x1 = x[1] - x[0] * x[0];
y1 = 1.0 - x[0];
y = 100.0 * x1 * x1 + y1 * y1;
return y;
}
// 関数３
double snx3(double *x)
{
double x1, y1, z1, y;
x1 = 1.5 - x[0] * (1.0 - x[1]);
y1 = 2.25 - x[0] * (1.0 - x[1] * x[1]);
z1 = 2.625 - x[0] * (1.0 - x[1] * x[1] * x[1]);
y = x1 * x1 + y1 * y1 + z1 * z1;
return y;
}