C言語によるファジィ(Fuzzy)推論コード

Go言語を使った自作のファジィ(Fuzzy)推論コードを探していたのですが、見つけることができずにちょっとショックを受けております。

というか、Go言語どころか、C言語でも見つけられなかった(いや、あるはずだ)ので、今、PCの中探して見つかったものをアップしておきます。

毎回スクラッチから作るのも迂遠なので。

↓のコードですが、こちらのコラムで使ったものだと思います。

/*
  このプログラムをfuzzy.cpp"という名前で保存して、
  gcc -g fuzzy.cpp -o fuzzy
  でコンパイルして下さい。
*/



#include <stdio.h>
#include <stdlib.h>


// 共通関数
double max(double a){
    return a;
}

double max(double a, double b){
    if (a > b) 
        return a;
    else 
        return b;
};

double max(double a, double b, double c ){
    return max(max(a,b), max(b,c)); 
};

double max(double a, double b, double c, double d ){
    return max(max(a,b,c), d); 
};

double min(double a){
    return a;
}

double min(double a, double b){
    if (b > a) 
        return a;
    else 
        return b;
};

double min(double a, double b, double c ){
    return min(min(a,b), min(b,c)); 
};

double min(double a, double b, double c, double d ){
    return min(min(a,b,c), d); 
};


// ファジィ表現
typedef enum scale {LESSLESS, LESS, ZERO, MORE, MOREMORE} SCALE;

// 前件部メンバーシップ関数(山3つ)クラス
class condition_MF3
{
private:
    double center;
    double width;
    SCALE express;
   
public:
    condition_MF3(double _center, double _witdth, SCALE _express){
        center = _center;
        width = _witdth;
        express = _express;

        // 使用できないファジィ表現を使った場合は止める        
        if ((express == LESSLESS) || (express == MOREMORE)){
            printf("wrong expression used \n");
            exit(0);
        }

    };
    double func(double _x);
};

double condition_MF3::func(double _x)
{
// x,yは、メンバーシップ関数上の座標を示す
    double x = _x;
    double y = 0.0; // yの値は、必ず0以上1以下になる

    if (express == LESS){
        if (x <= center - width){
            y = 1.0;
        }
        else if (x <= center){
            y = - 1.0 / width * (x - center);
        }
        else{
            y = 0.0;
        }
    }
    else if (express == ZERO){
        if (x <= center - width){
            y = 0.0;
        }
        else if (x <= center){
            y = 1.0 / width * (x - center) + 1.0;
        }
        else if (x <= center + width){
            y = -1.0 / width * (x - center) + 1.0;
        }
        else{
            y = 0.0;
        }
    }
    else if (express == MORE){
        if (x <= center){
            y = 0.0;
        }
        else if (x <= center + width){
            y = 1.0 / width * (x - center);
        }
        else{
            y = 1.0;
        }
    }
    else {
        printf("wrong expression\n");
        exit(1);
    }

    return y;
};

// 前件部メンバーシップ関数(山5つ)クラス
class condition_MF5
{
private:
    double center;
    double width;
    SCALE express;
   
public:
    condition_MF5(double _center, double _witdth, SCALE _express){
        center = _center;
        width = _witdth;
        express = _express;
    };
    double func(double _x);
};


double condition_MF5::func(double _x)
{
// x,yは、メンバーシップ関数上の座標を示す
    double x = _x;
    double y = 0.0; // yの値は、必ず0以上1以下になる

    if (express == LESSLESS){
        if (x <= center - 2.0 * width){
            y = 1.0;
        }
        else if (x <= center - width){
            y = - 1.0 / width * (x - (center - 2.0 * width)) + 1.0;
        }
        else{
            y = 0.0;
        }
    }
    else if (express == LESS){
        if (x <= center - 2.0 * width){
            y = 0.0;
        }
        else if (x <= center - width){
            y = 1.0 / width * (x - (center - width)) + 1.0;
        }
        else if (x <= center){
            y = -1.0 / width * (x - (center - width)) + 1.0; 
        }
        else{
            y = 0.0;
        }
    }
    else if (express == ZERO){
        if (x <= center - width){
            y = 0.0;
        }
        else if (x <= center){
            y = 1.0 / width * (x - center) + 1.0;
        }
        else if (x <= center + width){
            y = -1.0 / width * (x - center) + 1.0;
        }
        else{
            y = 0.0;
        }
    }
    else if (express == MORE){
        if (x <= center){
            y = 0.0;
        }
        else if (x <= center + width){
            y = 1.0 / width * (x - (center + width)) + 1.0;
        }
        else if (x <= center + 2.0 * width){
            y = -1.0 / width * (x - (center + width)) + 1.0; 
        }
        else{
            y = 0.0;
        }
    }
    else if (express == MOREMORE){
        if (x <= center + width){
            y = 0.0;
        }
        else if (x <= center + 2.0 * width){
            y = 1.0 / width * (x - (center + 2.0 * width)) + 1.0;
        }
        else{
            y = 1.0;
        }
    }

    return y;
};

// 後件部メンバーシップ関数(山3つ)クラス  
class action_MF3
{
private:
    double center;
    double width;
    SCALE express;

    double x;
    double y;

public:
    action_MF3(double _center, double _witdth, SCALE _express){

        y = 0.0; // yの値は、必ず0以上1以下になる
       
        center = _center;
        width = _witdth;
        express = _express;
        
        if (express == LESS){
            x = center - width;
        }
        else if (express == ZERO){
            x = center;
        }
        else if (express == MORE){
            x = center + width;
        }
        else{
            printf("wrong scale expression\n");
            exit(0);
        }
    };
    
    void func_Max(double b){
        y = max(b, y);
    };
    
// X座標を返す
    double func_X(void){
        return x;
    };
    
    // (最大値で更新された、最後の)Y座標を返す
    double func_Y(){
        return y;
    };

};

// 後件部メンバーシップ関数(山5つ)クラス  
class action_MF5
{
private:
    double center;
    double width;
    SCALE express;

    double x;
    double y;

public:
    action_MF5(double _center, double _witdth, SCALE _express){
        y = 0.0; // yの値は、必ず0以上1以下になる
       
        center = _center;
        width = _witdth;
        express = _express;
        
        if (express == LESSLESS){
            x = center - 2.0 * width;
        }
        else if (express == LESS){
            x = center - width;
        }
        else if (express == ZERO){
            x = center;
        }
        else if (express == MORE){
            x = center + width;
        }
        else if (express == MOREMORE){
            x = center + 2.0 * width;
        }
        else{
            printf("wrong scale expression\n");
        }
    };
    
    void func_Max(double b){
        y = max(b, y);
    };
    
// X座標を返す
    double func_X(void){
        return x;
    };
    
    // (最大値で更新された、最後の)Y座標を返す
    double func_Y(){
        return y;
    };

};

int main()
{
    // ダイエット経過期間
    condition_MF3 Period_Short(14, 14, LESS); // 0日
    condition_MF3 Period_Middle(14, 14, ZERO); // 14日
    condition_MF3 Period_Long(14, 14, MORE); // 28日

    // いわゆる「停滞期」(14日間あたりの平均減重変化量)
    condition_MF3 WeighChange_Small(1.0, 1.0, LESS); // 0kg
    condition_MF3 WeighChange_Middle(1.0, 1.0, ZERO); // 1kg
    condition_MF3 WeighChange_High(1.0, 1.0, MORE); // 2kg

    // 苦痛度
    condition_MF3 Suffer_Zerol(400.0, 400.0, LESS); // 0
    condition_MF3 Suffer_Middle(400.0, 400.0, ZERO); // 400
    condition_MF3 Suffer_Hight(400.0, 400.0, MORE); // 800

    // BMI
    condition_MF5 BMI_lowlow(22.0, 3.0, LESSLESS); // 痩せすぎ(15.0)
    condition_MF5 BMI_low(22.0, 3.0, LESS);// 痩せぎみ(19.0)
    condition_MF5 BMI_Normal(22.0, 3.0, ZERO);// 普通(22.0)
    condition_MF5 BMI_High(22.0, 3.0, MORE); // 太りぎみ(25.0)
    condition_MF5 BMI_HighHigh(22.0, 3.0, MOREMORE);// 太りすぎ(28.0)

    // 一日の摂取カロリー
    action_MF5 Eat_LittleLittle(2500, 1000, LESSLESS); // 500 Kcal
    action_MF5 Eat_Little(2500, 1000, LESS);// 1500 Kcal
    action_MF5 Eat_Normal(2500, 1000, ZERO);// 2500 Kcal
    action_MF5 Eat_Lot(2500, 1000, MORE);// 3500 Kcal
    action_MF5 Eat_LotLot(2500, 1000, MOREMORE);// 4500 Kcal



    
    
    // 超シンプルシミュレータ投入
    // 体重変動シミュレーションプログラム
    //基礎代謝量と求め方	?ハリス・ベネディクト方程式(日本人版) 計算式
    //http://www.kintore.info/kisotaisya_mass/index.cgi
    //【計算式】
    // 男性 66.5+(体重kg×13.8)+(身長cm×5.0)-(年齢×6.8)
    // 女性 665.1+(体重kg×9.6)+(身長cm×1.9)-(年齢×7.0)							
    // 江端智一 2014日2月4日現在のデータを使って計算
    double weight = 78.0;  // 開始時の体重 78.0kg
    double lengthw = 172.0;  // 身長 172cm
    double age = 49.0 + 94/365;  // 開始時の年齢 49歳と94日 (2014年2月4日現在)

    // 最初の基礎代謝カロリー 男性  66.5+(体重kg×13.8)+(身長cm×5.0)-(年齢×6.8)
    double consumption_calorie =  
        66.5 + weight * 13.8 + 172.0 * 5.0 - age * 6.8;
    // このカロリーに運動消費カロリーを加えて、
    // トータルの一日の消費カロリーを算出する

    // (現在の体重 66.5kgにおいて)一日運動消費カロリー 708kcalと推定し、
    // このカロリーは、体重に比例するものと仮定する。
    consumption_calorie += weight/66.5 * 708.0;
    printf("consumption_calorie = %f\n",  consumption_calorie);

    double last_weight[14]; // 過去14日分の体重
    
    // 最初は全部同じ体重とする 14日経過までは同じ体重とする
    for (int i = 0; i++; i < 14){
        last_weight[i] = weight;
    }

    // ここからループ開始
    for (int day=0; day < 1460; day++){ // 1460日は、ちょうど4年分の日数
        // 7kcal = 1g = 0.001kg とする    

        // 体重データを一日分ずらす
        for (int i = 0; i++; i < 14-1 ){
            last_weight[i+1] = last_weight[i];
        }
        last_weight[0] = weight;
        
        // 線形補完の傾きだけを求める(14日分)
        // http://d.hatena.ne.jp/rainlib/20090112/1231735459
        double k1=0,k2=0,k3=0,k4=0;
        for (int i = 0; i++; i < 14 ){
            k1 += (double)i * (double)i;
            k2 += (double)i * last_weight[i];
            k3 += (double)i;
            k4 += last_weight[i] *last_weight[i];
        }
        double dW = (14.0 * k2 - k3 * k4) / (14.0 * k4 - k3*k3); // 14

        // BMIの計算式 体重(Kg)÷身長(m)÷身長(m) 江端の身長172cm
        double bmi = weight / 1.72 / 1.72;
        
        // [ルール1] ダイエット開始直後は、「無条件にがんばれる」
        double r1 = min(Period_Short.func(day));
        Eat_Little.func_Max(r1);
        //printf("r1: %f\n",r1);
        
        //printf("Eat_Little.func_X() = %f\n",Eat_Little.func_X());
        //printf("Eat_Little.func_Y() = %f\n",Eat_Little.func_Y());
        
        // [ルール2] ダイエット一ヶ月を過ぎても、停滞し続けると「切れる」
        double r2 = min(Period_Long.func(day),WeighChange_Small.func(dW));
        Eat_LotLot.func_Max(r2);
        //printf("r2: %f\n",r2);
        
        // [ルール3] ダイエット2週間を過ぎて、停滞し続けると「緩んでくる」
        double r3 = min(Period_Middle.func(day),WeighChange_Small.func(dW));
        Eat_Normal.func_Max(r3);
        //printf("r3: %f\n",r3);
        
        // [ルール4] 停滞期がなく順調に体重が落ちている場合は「がんばれる」
        double r4 = min(WeighChange_High.func(dW));
        Eat_Little.func_Max(r4);
		//        printf("r4: %f\n",r4);
        
        // [ルール5] ダイエット一ヶ月を過ぎて、再び太り出してくると、またダイエットを再開してしまう
        double r5 = min(Period_Long.func(day),BMI_High.func(bmi));
        Eat_Little.func_Max(r5);
        //printf("r5: %f\n",r5);
        
        // [ルール6] ダイエット一ヶ月を過ぎて、太り過ぎると、滅茶苦茶なダイエットを始めてしまう
        double r6 = min(Period_Long.func(day),BMI_HighHigh.func(bmi));
        Eat_LittleLittle.func_Max(r6);

        //printf("r6: %f\n",r6);
		/*        
        printf("Eat_LittleLittle.func_X() = %f\n",Eat_LittleLittle.func_X());
        printf("Eat_LittleLittle.func_Y() = %f\n",Eat_LittleLittle.func_Y());
        printf("Eat_Little.func_X() = %f\n",Eat_Little.func_X());
        printf("Eat_Little.func_Y() = %f\n",Eat_Little.func_Y());
        printf("Eat_Normal.func_X() = %f\n",Eat_Normal.func_X());
        printf("Eat_Normal.func_Y() = %f\n",Eat_Normal.func_Y());
        printf("Eat_Lot.func_X() = %f\n",Eat_Lot.func_X());
        printf("Eat_Lot.func_Y() = %f\n",Eat_Lot.func_Y());
        printf("Eat_LotLot.func_X() = %f\n",Eat_LotLot.func_X());
        printf("Eat_LotLot.func_Y() = %f\n",Eat_LotLot.func_Y());
        */
        
        double source_calorie =  
            (Eat_LittleLittle.func_X() * Eat_LittleLittle.func_Y() + 
             Eat_Little.func_X() * Eat_Little.func_Y() +
             Eat_Normal.func_X() * Eat_Normal.func_Y() +
             Eat_Lot.func_X() * Eat_Lot.func_Y() +
             Eat_LotLot.func_X() * Eat_LotLot.func_Y() )  
            /
            (Eat_LittleLittle.func_Y() + 
             Eat_Little.func_Y() +
             Eat_Normal.func_Y() +
             Eat_Lot.func_Y() +
             Eat_LotLot.func_Y() ) ;
        
        printf("\n source_calorie = %f\n",source_calorie);

        double diff_weight = (source_calorie - consumption_calorie)/7.0 / 1000.0;
        weight += diff_weight;
        age += 1.0/365.0;
        consumption_calorie =  66.5 + weight * 13.8 + 172.0 * 5.0 - age * 6.8;
        consumption_calorie += weight/66.5 * 708.0;
        //       printf("day:%d\tweight = %f\tconsumption_calorie =%f\n", day,weight, consumption_calorie);
    }



}

ファジィ推論

 

2022/09,江端さんの技術メモ

Posted by ebata