Matlab使用m文件调用simulink（附matlab代码）

B站讲解视频
我就是按这个视频写的代码
我用simulink打了一个简单地模型，希望用遗传算法识别模型中的参数。
实现这个目的得流程是将simulink模型中的参数设为变量，每运行一次模型就将结果和对应的参数值读取到matlab中，然后通过遗传算法计算更优参数值，再将参数值赋给模型重新计算，得到结果开始循环。
假如不想自己搭模型可以直接用我的，
链接：https://pan.baidu.com/s/1di1_5vJIZXUdWTOsw_EPvA
提取码：7788
–来自百度网盘超级会员V5的分享
如果要自己搭模型的话记得按视频里的操作进行设置。

我真的好心烦csdn传资源这个事，这不是一个分享平台吗，我把我的学习资源传上来供大家学习，结果其他人想下载还要充会员，甚至有的充了会员还要花钱下载。有时候看一篇博客，博主说欢迎大家下载他的文件，大家一起学习交流，结果一点开需要花钱才能下，博主自己知道自己的资源需要别人花钱才能下载吗，给博主分成是多少啊。
要是博主自己要收钱那无可厚非，可我的资源是想无偿给大家的啊，而且我一个资源分享者都无法评论我的资源，要下载后才能评论，真是离谱，一个好端端的分享平台竟然干这事。
所以我把资源放在百度网盘里了，大家可以免费下载。

评论说链接挂掉了，很奇怪我明明选的是永久有效，我再放一个新的
链接：https://pan.baidu.com/s/1Bx5lbmHrHxEBMFJTySp7TQ
提取码：qtyt
–来自百度网盘超级会员V6的分享

然后这是使用的matlab代码

clc
close all
clear
mdlName = 'model';
%参数标准值
b1_ref = 1;
b2_ref = 1;
a11_ref = 0.4;
a12_ref = 0.4;
a21_ref = 0.2;
a22_ref = 0.2;

b1_test = 0;
b2_test = 0;
a11_test = 0;
a12_test = 0;
a21_test = 0;
a22_test = 0;


% 参数范围,按照b1,b2,a11,a12,a21,a22的顺序
lb = [0,0,0,0,0,0];
ub = [2,2,1,1,0.5,0.5];

% 种群数量
popSize = 50;

% 迭代次数
maxGen = 100;

% 交叉率
crossRate = 0.9;

% 变异率
mutationRate = 0.1;

% 初始化种群
pop = rand(popSize,6).*(ub-lb)+lb;
fitness = zeros(popSize,maxGen+1);
n = 0;%记录最优适应度相同的次数
fitnessmax = 0;%记录最优适应度
% 迭代
for i = 1:maxGen
    mutationRate = 0.1;
    
    
    for j = 1:popSize
        b1 = pop(j,1);
        b2 = pop(j,2);
        a11 = pop(j,3);
        a12 = pop(j,4);
        a21 = pop(j,5);
        a22 = pop(j,6);

        % 评估适应度
        load_system(mdlName);
        cs = getActiveConfigSet(mdlName);
        model_cs = cs.copy;
        simOut = sim(mdlName,model_cs);
        x1 = simOut.logsout{1}.Values.Data;
        x2 = simOut.logsout{2}.Values.Data;
        x1_ref = simOut.logsout{3}.Values.Data;
        x2_ref = simOut.logsout{4}.Values.Data;
        fitness(j,i) = fun(x1,x1_ref,x2,x2_ref);
    end
    
    % 选择操作
    [~,index] = sort(fitness(:,i));
    if fitnessmax == fitness(index(1))
        n = n+1;
    else
        n = 0;
    end
    fitnessmax = fitness(index(1));%更新最佳适应度
    parents = pop(index(1:popSize/2),:);
    
    % 交叉操作
    offsprings = zeros(popSize/2,6);
    for j = 1:popSize/2
        if rand < crossRate
            parent1 = parents(randi(popSize/2),:);
            parent2 = parents(randi(popSize/2),:);
            nums = randperm(6, 3);
            %将两个亲代的三个随机基因进行交换
            for m = 1:length(nums)
                parent1(nums(m)) = parent2(nums(m));
            end
            offsprings(j,:) = parent1;
        else
            offsprings(j,:) = parents(randi(popSize/2),:);
        end
    end
    
    % 变异操作
    if fitnessmax >20 && n == 7%n=7设为阈值
        mutationRate = 0.9;%增大变异概率，进入物种大爆发阶段
        n = 0;
        for k = 5:popSize/2-5
            c = ceil(rand(1,1)*6);
            parents(k,c) = rand(1,1).*(ub(c)-lb(c))+lb(c);
        end
        for k = popSize/2-5:popSize/2
            parents(k,:) = rand(1,6).*(ub-lb)+lb;
        end
    end
    for j = 1:popSize/2
        if rand < mutationRate
            c = ceil(rand(1,1)*6);
            offsprings(j,c) = rand(1,1).*(ub(c)-lb(c))+lb(c);
        end
    end
    offsprings(end,:) = rand(1,6).*(ub-lb)+lb;
    % 合并新一代种群
    pop = [parents;offsprings]; 
%     fitness(:,i)
     i
%     n
end

% 返回最优解

 for j = 1:popSize
        b1 = pop(j,1);
        b2 = pop(j,2);
        a11 = pop(j,3);
        a12 = pop(j,4);
        a21 = pop(j,5);
        a22 = pop(j,6);

        % 评估适应度
        load_system(mdlName);
        cs = getActiveConfigSet(mdlName);
        model_cs = cs.copy;
        simOut = sim(mdlName,model_cs);
        x1 = simOut.logsout{1}.Values.Data;
        x2 = simOut.logsout{2}.Values.Data;
        x1_ref = simOut.logsout{3}.Values.Data;
        x2_ref = simOut.logsout{4}.Values.Data;
        fitness(j,end) = fun(x1,x1_ref,x2,x2_ref);
end

[~,index] = min(fitness(:,end));
x = pop(index,:);
b1_test = x(1);
b2_test = x(2);
a11_test = x(3);
a12_test = x(4);
a21_test = x(5);
a22_test = x(6);
load_system(mdlName);
cs = getActiveConfigSet(mdlName);
model_cs = cs.copy;
simOut = sim(mdlName,model_cs);


plot(fitness(1,:),'LineWidth',3)
xlabel('迭代次数')
ylabel('适应度值')

不好意思忘记放里面适应度的计算函数了，其实就是误差啦。这是fun函数的写法。

function e =fun(x1,x1_ref,x2,x2_ref)
    e1 = (sum((x1-x1_ref).*(x1-x1_ref)))^(0.5);
    e2 = (sum((x2-x2_ref).*(x2-x2_ref)))^(0.5);
    e = (e1+e2);
end