Trans论文复现：含碳交易的动态直流最优潮流ADMM优化方法程序代码！

本程序参考期刊论文《A Distributed Dual Consensus ADMM Based on Partition for DC-DOPF with Carbon Emission Trading》，本程序算例丰富，注释清晰，质量很高，干货满满，小编非常推荐，下面对论文和程序作简单介绍。

本文提出了含碳排放交易的直流动态最优潮流（DC-DOPF-CET）分布式交替方向（ADMM）优化方法。常见的ADMM方法公开相邻子系统之间的边界总线和分支信息，本文所提出的方法（DC-ADMM-P）采用了一种新的策略，该策略使用一致性ADMM来解决DC-DOPF-CET的对偶问题，同时只公开相邻子系统之间的边界分支信息。此外，通过减少双乘法的数量和采用改进的更新步骤，DC-ADMM-P的收敛性有很大改善。DC-ADMM-P在6-1062系统上进行了测试，并与其他分布式/去中心化方法进行了比较。仿真结果验证了DC-ADMM-P在求解复杂（非线性）的DC-DOPF问题时的高效性，该问题可以被公式化为凸可分离函数。同时，它还通过减少对偶乘法的数量和采用新的乘法更新策略来提高收敛性。

程序组成：

程序由多个子函数组成，有几个子函数可以单独运行，计算算例中的成本、碳交易成本、算法时间。

程序结果：

本图片为程序在matlab上的运行结果，在运行程序时，求解次数与相应的成本对应，可以更好的看出求解结果变化情况。

以上图片与论文中的相对应，展示复现精度。

部分程序：

%% ***************************************************************
function [  ] = DCOPF_ADMM(  )
            %---------------------------------Console 控制台-------------------------------------------start-------------------------------%
%             FileName = 'SCUC_dat/SCUC6.txt';    %Corresponding to the 6-bus System;对应文中6bus例子
%             FileName = 'SCUC_dat/SCUC30.txt'; %Corresponding to the 30-bus System;对应文中30bus例子
%             FileName = 'SCUC_dat/SCUC6-2.txt';  %Corresponding to the 6-bus System in literature [7];对应文献[7]中的6-bus system
%             FileName = 'SCUC_dat/SCUC1062-2.txt';  %Corresponding to the 1062-bus System;对应文中1062-bus例子
%             FileName = 'SCUC_dat/RTS48.txt'; %Corresponding to the RTS-48 bus system.The test system can obtain from [44];对应文中RTS0-48 bus例子
            FileName = 'SCUC_dat/DDOPF118.txt'; %Corresponding to the 118-bus system in literature [7]；对应文献[7]中的118-bus system
            type_of_pf = 'DC';  %The type of optimal power flow. 
            useImprovedPratition = 'yes'; %Reducing the number of dual multipliers. Choose 'yes' or 'no'. Corresponding to the Section III-C.          
            useImprovedMultiplier = 'no';  %Improving the update step of the multipliers. Choose 'yes' or 'no'. Corresponding to the Section III-E.
            includeDR = 'yes'; %Considering demand response. Choose 'yes' or 'no'. Corresponding to the Section IV-D.
            useParfor = 'no';  %Using parfor in Matlab. Choose 'yes' or 'no'.
            includeCet = 'no'; %Considering carbon emission trading. Choose 'yes' or 'no'.
            termination_criteria_type = 1; %The termination criteria strategy for ADMM. Choose '1' or '2'. Corresponding to the S_C 1 and S_C 2 in the second paragraph of Section IV.
            partionMethod = 3; %The partition strategy. For example, the Subsystems An in Table VI are the partions of 30-bus system. The implementation code is in the function partitionNode( ... ). 
            parallel_worker = 4; %The number of worker  (parfor in Matlab).
            M = 4000; %Maximum number of iterations for ADMM.
            rou = 0.31; %Penalty parameter
            TOL = 1e-5; %Tolerance parameter
            xishu = 0.54; %Undertermined relaxation factor in (38).
            
            %Corresponding to the Section IV-D.
            if isequal(includeDR,'yes')
                parameter.priceElasticity = -1; %Price elasticity [48].
                parameter.piecewiseNumber = 8; %An eight-piece piece wise segments is utilized to approximate the preice-demand curves.
                parameter.referencePoint = [80,25];%The reference point for the demand response curve is ( reference demand, reference price);需求曲线的参考点坐标，[需求负荷,参考价格]
            end

            %Corresponding to S_C 2 in the manuscript. The definition of S_C 2 is obtained from [7].
            if isequal(termination_criteria_type,2) 
                rk_TOL = 0.01; %Tolerance of primal residual,default value is 0.01 
                sk_TOL = 0.01; %Tolerance of dual residual,default value is 0.01 
            end

            %---------------------------------Console 控制台-------------------------------------------end-------------------------------%

            %Read datas
            FileIndex = strfind(FileName,'/');
            if isequal(FileName(FileIndex+1:FileIndex+5),'DDOPF') 
                %Corresponding to the 'SCUC_dat/DDOPF.......txt' files;DDOPF数据不考虑CET，但考虑DR；；不存在一个节点上有多个机组的情况
                isRTS = 'no'; %默认值,default value is 'no'
                includeDR = 'yes'; %默认值,default value is 'yes'
                includeCet = 'no'; %默认值,default value is 'no'
                parameter.includeDR = includeDR;
                parameter.includeCet = includeCet;
                parameter.isRTS = isRTS;
                SCUC_data = ReadDataDCDOPF(FileName,parameter);
            elseif isequal(FileName(FileIndex+1:FileIndex+3),'RTS')
                %%Corresponding to the 'SCUC_dat/RTS.......txt' files;RTS数据不考虑CET和DR；存在一个节点上有多个机组的情况
                isRTS = 'yes';  %是否使用RTS数据，'yes'为使用RTS数据,'no'为不适用RTS数据。,default value is 'yes'
                includeDR = 'no'; %默认值,default value is 'no'
                includeCet = 'no'; %默认值,default value is 'no'
                parameter.isRTS = isRTS;
                parameter.includeDR = includeDR;
                parameter.includeCet = includeCet;
                SCUC_data = ReadDataDCDOPF(FileName,parameter);
            elseif isequal(FileName(FileIndex+1:FileIndex+4),'SCUC')%SCUC数据
                %%Corresponding to the 'SCUC_dat/SCUC.......txt' files;本文SCUC数据，不考虑DR；不存在一个节点上有多个机组的情况
                isRTS = 'no'; %默认值,default value is 'no'
                includeDR = 'no'; %默认值,default value is 'no'
                parameter.includeDR = includeDR;
                parameter.includeCet = includeCet;
                SCUC_data = ReadDataSCUC(FileName,parameter);
            end
            
            T = SCUC_data.totalLoad.T;  % The number of time periods. 时段数T
            G = SCUC_data.units.N;      % The total number of units.发电机数
            N = SCUC_data.baseparameters.busN;  % The total number of buses.节点总数

            %The partition in test system; 下面的代码生成{1,...,N}的一个划分 PI
            best_partition = 1;     % Default value is 1;=1 直接按照最好的划分方式划分   =0，则按照给定n来划分
            [PI, allNode,PINumber] = partitionNode(SCUC_data.baseparameters.busN,N,best_partition,partionMethod);
            n = size(PI,1);  %The number of partition set; 划分的子集个数

            %Calculate the network admittance matrix;获取潮流节点导纳矩阵B
            Y = SCUC_nodeY(SCUC_data, type_of_pf);
            B = -Y.B;
            

从以上程序部分可以看出，控制台中可以手动调节想要的算例，简单易操作。展示到这里，您看出该代码的精美之处了吗，欢迎关注并获取完整版程序代码哦~