moead算法c++代码

A：以下是MOEA/D算法的C代码。

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <float.h>
#include <string.h>

#define N 1000 //population size
#define D 7 //number of decision variables
#define M 2 //number of objectives
#define NR 1000 //number of iterations
#define T 20 //neighbourhood size
#define CR 1.0 //crossover probability
#define F 0.5 //differential weight
#define EPS 1e-6 //epsilon

double genrand(){ //random number generator; Mersenne Twister from Matsumoto and Nishimura (1998)
    unsigned long x[624], y; 
    static int p=0; 
    static unsigned long mag01[2]={0x0UL,0x9908b0dfUL}; 
    int i; 

    if (p == 0) 
    for (i=0;i<624;i++) { 
        y = (x[i] & 0x80000000UL) + (x[(i+1)%624] & 0x7fffffffUL); 
        x[i] = x[(i + 397) % 624] ^ (y >> 1) ^ mag01[y & 0x1UL]; 
    } 

    y = x[p]; 
    p = (p + 1) % 624; 

    y ^= (y >> 11); 
    y ^= (y << 7) & 0x9d2c5680UL; 
    y ^= (y << 15) & 0xefc60000UL; 
    y ^= (y >> 18); 

    return (double)y/(double)((unsigned long)0xffffffffUL+1.0); 
} 

void set_ref(){ //initialize reference point
    int i;

    for(i=0; i<M; i++) f_ref[i] = DBL_MAX;

    i=0;
    while(i<N){
        double x[D];
        int j;
        for(j=0; j<D; j++) x[j] = 2.0*genrand() - 1.0;
        decode(x, N,P,Q,nvar, node, edge); //decode
        evaluate(x, &f[i*M]); //evaluate
        int flag=0;
        int k;
        for(k=0;k<i;k++){
            int l=0;
            while(l<M && f_ref[l] - f[k*M+l]>EPS) l++;
            if(l==M){
                flag=1;
                break;
            }
        }
        if(flag==0){
            for(k=0;k<M;k++) f_ref[k]=f[k];
            i++;
        }
    }
}

double sgn(double x){
    if(x<0.0) return -1.0;
    if(x>0.0) return 1.0;
    return 0.0; 
}

double power(double x, int n){
    double y=1.0;
    int i;
    for(i=0;i<n;i++)
        y *= x;
    return y;
}

//----------------------------------------------------------------------------

void func(vector x, vector f){ //multi-objective functions
    //f1 = x1
    //f2 = g(x)h(f1,g(x))
    double a=1.0;
    double b=1.0;

    int i;
    for(i=2;i<D;i++){
        a += x[i]*x[i] - 10.0*cos(4.0*PI*x[i]);
    }

    f[0] = x[0];
    f[1] = a*(1.0-sqrt(x[0]/a)-x[0]*sin(10.0*PI*x[0]))*b;
}

void evaluate(double x[], double f[]){ //evaluate an individual
    vector vx = alloc_vector(D);
    vector vf = alloc_vector(M);
    int i;
    for(i=0; i<D; i++) vx[i] = x[i];
    func(vx,vf);
    for(i=0; i<M; i++) f[i] = vf[i];
    free_vector(vf);
    free_vector(vx);    
}

//----------------------------------------------------------------------------

void crossover(double xp[], double x1[], double x2[], double x3[]){ //crossover
    int idx = (int)(genrand()*D);
    int i;
    for(i=0; i<D; i++){
        if(genrand()<CR || i == idx){
            x1[i] = xp[i] + F*(x2[i] - x3[i]);
            if(x1[i]<-1.0) x1[i] = -1.0;
            if(x1[i]>1.0) x1[i] = 1.0;
        }else{
            x1[i] = xp[i];
        }
    }
}

void mutation(double xp[], double x1[]){ //mutation
    int i;
    for(i=0; i<D; i++){
        x1[i] = xp[i];
        if(genrand()<0.1){
            x1[i] = 2.0*genrand() - 1.0;
        }
    }
}

void select_parent(int p[], int ci){ //select a parent based on binary tournament selection
    int i;
    for(i=0; i<2; i++){
        int j = (int)(genrand()*N);
        if(f[j*M+ci]<f[p[i]*M+ci]){
            p[i] = j;
        }
    }
}

void moead(){ //MOEA/D -- main procedure
    int t, i, j, k;
    for(i=0;i<N;i++)
        for(j=0;j<M;j++)
            f[i*M+j]=0.0;
    set_ref(); //initialize reference point

    for(t=0; t<NR; t++){ 
        printf("# %d\n",t);

        for(i=0; i<N; i++){ //mutate every individual
            int p[2] = {i, i};
            select_parent(p, 0);
            crossover(x[p[0]], xp, x[p[1]], x[perm[i][genrand()*T]]);
            mutation(xp, x1);
            decode(x1, N,P,Q,nvar, node, edge);
            evaluate(x1, f1);

            for(j=0; j<T; j++){ //update neighbourhood individuals
                k = nei[i][j];
                double nf[M];
                evaluate(x[k], nf);
                int flag=0;
                int l=0;
                while(l<M && nf[l]-f[k*M+l]<EPS) l++;
                if(l==M){
                    double dis = 0.0;
                    int p;
                    for(p=0;p<D;p++) dis += pow(x1[p]-x[k][p],2.0);
                    if(dis < eps){
                        memcpy(x[k], x1, sizeof(double)*D);
                        memcpy(f+k*M, f1, sizeof(double)*M);
                    }else if(genrand()<0.5){
                        memcpy(x[k], x1, sizeof(double)*D);
                        memcpy(f+k*M, f1, sizeof(double)*M);
                    }
                }
            }
        }
    }
}

//----------------------------------------------------------------------------

int main(int argc, char *argv[]){
    init();
    moead();
    return 0;
}

请注意这是一个示例代码，并且可能不完全符合您的实际需求。