请修改这段代码，使得运行过程中指数函数的指数部分只有常数
# -*- coding: UTF-8 -*-
import os
import pandas as pd
import math
import random
import operator
from deap import creator, base, tools, gp, algorithms
import numpy
from fontTools.misc.py23 import xrange
from matplotlib import pyplot as plt
import pygraphviz as pgv
import datetime
import sys
path = os.getcwd()
curr_time = datetime.datetime.now()
time_str = datetime.datetime.strftime(curr_time, '%Y%m%d%H%M')
def text_save(filename, data):
file = open(filename, 'a')
for i in range(len(data)):
s = str(data[i]).replace('[', '').replace(']', '').replace('\n', '')
s = s.replace("'", '').replace(',', '') + '\n'
file.write(s)
file.close()
print("Saving the file succeeded")
return
def ProDiv(left, right):
if right == 0:
return 1000000
else:
return left / right
def ProSqrt(x):
if x <= 0:
return 1000000
else:
return x ** 0.5
def Prolog(x):
if x <= 0:
return 1000000
else:
return math.log(x)
def Proexp(x):
if x >= 100:
return 1000000
else:
return math.exp(x)
def Prosquare(x):
if x <= 0:
return 1000000
else:
if x ** 2 > 1e30:
return 1000000
else:
return x ** 2
def Propow(down, up):
"""
安全的指数计算函数，指数必须是纯常数且范围在(-10,10)
参数:
down: 底数（可以包含变量）
up: 指数（必须是纯常数，范围在-10到10之间）
返回:
计算结果，错误时返回1000000
"""
# 确保指数是有效常数
# 1. 检查是否为数字类型
if not isinstance(up, (int, float)):
return 1000000
# 2. 检查是否为特殊值（NaN或无穷大）
if math.isnan(up) or math.isinf(up):
return 1000000
# 3. 检查指数范围（-10 < up < 10）
if up <= -10 or up >= 10:
return 1000000
# 处理底数为零的情况
if down == 0:
return 0.0 if up > 0 else 1000000 # 0的正指数=0，负指数未定义
# 处理负底数
if down < 0:
return 1000000
# 使用对数预测结果大小，防止溢出
try:
# 计算结果的绝对值的对数
log_value = abs(up) * math.log10(abs(down))
# 检查结果是否过大(>10^50)或过小(<10^-50)
if log_value > 50 or log_value < -50:
return 1000000
except (ValueError, OverflowError):
return 1000000
# 正常计算并捕获可能的溢出错误
try:
result = down ** up
# 额外检查是否得到无穷大
if math.isinf(result):
return 1000000
return result
except (OverflowError, ValueError):
return 1000000
# 更可靠的常数节点检查
def is_constant_node(node):
"""检查节点是否是真正的常数节点"""
# 常数节点必须有 value 属性
if not hasattr(node, 'value'):
return False
# 确保值是可接受的数字类型
if not isinstance(node.value, (int, float)):
return False
# 确保值在合理范围内
if abs(node.value) > 1e10: # 防止过大或过小的常数
return False
return True
# 更严格的幂函数检查
def check_pow_arguments(individual):
"""确保所有 Propow 节点的指数参数都是常数"""
# 首先将个体转换为字符串表达式
expr_str = str(individual)
# 查找所有 Propow 函数调用
pow_calls = []
start_idx = 0
while True:
idx = expr_str.find("Propow(", start_idx)
if idx == -1:
break
# 找到匹配的右括号
depth = 1
end_idx = idx + len("Propow(")
while end_idx < len(expr_str) and depth > 0:
if expr_str[end_idx] == '(':
depth += 1
elif expr_str[end_idx] == ')':
depth -= 1
end_idx += 1
# 提取参数部分
args_str = expr_str[idx + len("Propow("):end_idx - 1].strip()
args = [arg.strip() for arg in args_str.split(',') if arg.strip()]
# 确保有两个参数
if len(args) == 2:
pow_calls.append((args[0], args[1]))
start_idx = end_idx
# 检查每个幂调用的指数参数
for base, exponent in pow_calls:
# 指数必须是数字（常数）
try:
float(exponent)
except ValueError:
# 如果无法转换为浮点数，说明是变量
print(f"发现非法幂函数表达式: Propow({base}, {exponent})")
return False
return True
# 在初始种群生成时强制验证
def create_valid_individual():
"""创建并验证个体，直到找到有效的个体"""
max_attempts = 100
attempts = 0
while attempts < max_attempts:
ind = tool.individual()
if check_pow_arguments(ind):
return ind
attempts += 1
# 如果多次尝试后仍无法生成合法个体，创建一个安全的默认个体
# 例如: Propow(1.0, 1.0)
default_expr = gp.PrimitiveTree.from_string(
"Propow(1.0, 1.0)", pset
)
return creator.Individual(default_expr)
# 修改种群创建方法
def create_valid_population(n):
"""创建并验证种群"""
pop = []
for _ in range(n):
pop.append(create_valid_individual())
return pop
# 增强的交叉操作
def custom_cxOnePoint(ind1, ind2):
max_attempts = 50
attempts = 0
while attempts < max_attempts:
new_ind1, new_ind2 = gp.cxOnePoint(ind1, ind2)
if check_pow_arguments(new_ind1) and check_pow_arguments(new_ind2):
return new_ind1, new_ind2
attempts += 1
# 如果多次尝试后仍无法生成合法个体，返回原始个体
return ind1, ind2
# 增强的变异操作
def custom_mutate(individual, expr, pset):
max_attempts = 50
attempts = 0
while attempts < max_attempts:
new_individual, = gp.mutUniform(individual, expr, pset)
if check_pow_arguments(new_individual):
return new_individual,
attempts += 1
# 如果多次尝试后仍无法生成合法个体，返回原始个体
return individual,
def Psin(l, r):
return math.sin(l * r)
def Pcos(l, r):
return math.cos(l * r)
def generate_random_float(min_value, max_value, decimal_places):
scale = 10 ** decimal_places
random_float = random.uniform(min_value, max_value)
rounded_float = round(random_float, decimal_places)
scaled_float = int(rounded_float * scale) / scale
return scaled_float
min_val = 0.0
max_val = 1.0
decimal_places = 4
# 定义 primitive set
pset = gp.PrimitiveSet('main', 4) # 定义变量个数
pset.renameArguments(ARG0='Fr') # 定义变量名称
pset.renameArguments(ARG1='We')
pset.renameArguments(ARG2='Re')
pset.renameArguments(ARG3='Bo')
# 只保留乘除和幂运算
pset.addPrimitive(operator.mul, 2, name='mul')
pset.addPrimitive(ProDiv, 2, name='ProDiv')
pset.addPrimitive(Propow, 2, name='Propow')
# 添加常数生成器
pset.addEphemeralConstant('c1', lambda: generate_random_float(-2, 2, 4))
pset.addEphemeralConstant('c2', lambda: generate_random_float(-1, 1, 4))
pset.addEphemeralConstant('c3', lambda: generate_random_float(0, 1, 4))
pset.addEphemeralConstant('c4', lambda: generate_random_float(1, 2, 4))
# 创建fitness类、individual类
creator.create('FitnessMin', base.Fitness, weights=(-1.0,))
creator.create('Individual', gp.PrimitiveTree, fitness=creator.FitnessMin)
# 定义个体的生成方法、种群的生成方法
tool = base.Toolbox()
tool.register('expr', gp.genHalfAndHalf, pset=pset, min_=1, max_=3)
tool.register('individual', tools.initIterate, creator.Individual, tool.expr)
tool.register('population', tools.initRepeat, list, tool.individual)
tool.register('compile', gp.compile, pset=pset)
tool.register('expr_mut', gp.genFull, pset=pset, min_=0, max_=2) # 生成一个subtree
# 修改变异操作，确保变异后的个体符合要求
def custom_mutate(individual, expr, pset):
new_individual, = gp.mutUniform(individual, expr, pset)
while not check_pow_arguments(new_individual):
new_individual, = gp.mutUniform(individual, expr, pset)
return new_individual,
# 修改交叉操作，确保交叉后的个体符合要求
def custom_cxOnePoint(ind1, ind2):
new_ind1, new_ind2 = gp.cxOnePoint(ind1, ind2)
while not check_pow_arguments(new_ind1) or not check_pow_arguments(new_ind2):
new_ind1, new_ind2 = gp.cxOnePoint(ind1, ind2)
return new_ind1, new_ind2
# 定义评价函数
def fit_evaluation(individual, input1, input2, input3, input4, input5, output):
func = tool.compile(expr=individual) # pset上面已经给过了
sqerrors = ((abs((1 + func(Fr, We, Re, Bo)) * dPa - dPe) / dPe * 100) for Fr, We, Re, Bo, dPa, dPe in
zip(input1, input2, input3, input4, input5, output))
return math.fsum(sqerrors) / len(output), # 必须返回一个tuple
path = os.getcwd()
filename = path + '/' + 'data_horizontal_mg_4g' + '.xlsx'
book = pd.read_excel(filename)
data = book.values
# 检查文件是否存在
if os.path.exists(filename): # 确保已经导入了os模块
try:
# 读取Excel文件
book = pd.read_excel(filename)
data = book.values
print(f"成功读取文件: {filename}")
print(f"数据形状: {data.shape}")
except Exception as e:
print(f"读取文件时发生错误: {e}")
else:
print(f"错误：文件不存在 - {filename}")
# 可以在这里添加创建文件或其他处理逻辑
# x = data[:, 37]
# X = (data[:, 31]/data[:, 32])**0.5
Fr = data[:, 28] # (data[:, 1]+2.2)
We = data[:, 32]
Re = data[:, 24]
Bo = data[:, 36]
# N=(WeG*Fr**2)**0.25
# Pred = data[:, 19]
# f = data[:, 35]/data[:, 36]
dPl = data[:, 61]
dPe = data[:, 65]
Y = (data[:, 67] / data[:, 66]) ** 0.5
dPa = data[:, 66] # /1000000000
# 定义evaluate、select、mate、mutate（这几个名字必须这样取，否则出错）
tool.register('evaluate', fit_evaluation, input1=Fr, input2=We, input3=Re, input4=Bo, input5=dPa, output=dPe)
tool.register('select', tools.selTournament, tournsize=3)
# 注册修改后的交叉操作
tool.register('mate', custom_cxOnePoint)
# 注册修改后的变异操作
tool.register('mutate', custom_mutate, expr=tool.expr_mut, pset=pset)
# 限制树深度
tool.decorate("mate", gp.staticLimit(key=operator.attrgetter("height"), max_value=10))
tool.decorate("mutate", gp.staticLimit(key=operator.attrgetter("height"), max_value=10))
def drawBestInd(expr): # 输入参数是一个表达式（expr/tree），这个例子中就是individual
nodes, edges, labels = gp.graph(expr)
g = pgv.AGraph()
g.add_nodes_from(nodes)
g.add_edges_from(edges)
g.layout(prog="dot")
for i in nodes:
n = g.get_node(i)
n.attr["label"] = labels[i]
# curr_time = datetime.datetime.now()
# time_str = datetime.datetime.strftime(curr_time, '%Y%m%d%H%M')
g.draw(time_str + "tree.pdf")
def varOr(population, toolbox, lambda_, cxpb, mutpb):
"""Part of an evolutionary algorithm applying only the variation part
(crossover, mutation **or** reproduction). The modified individuals have
their fitness invalidated. The individuals are cloned so returned
population is independent of the input population.
:param population: A list of individuals to vary.
:param toolbox: A :class:`~deap.base.Toolbox` that contains the evolution
operators.
:param lambda\_: The number of children to produce
:param cxpb: The probability of mating two individuals.
:param mutpb: The probability of mutating an individual.
:returns: The final population.
The variation goes as follow. On each of the *lambda_* iteration, it
selects one of the three operations; crossover, mutation or reproduction.
In the case of a crossover, two individuals are selected at random from
the parental population :math:`P_\mathrm{p}`, those individuals are cloned
using the :meth:`toolbox.clone` method and then mated using the
:meth:`toolbox.mate` method. Only the first child is appended to the
offspring population :math:`P_\mathrm{o}`, the second child is discarded.
In the case of a mutation, one individual is selected at random from
:math:`P_\mathrm{p}`, it is cloned and then mutated using using the
:meth:`toolbox.mutate` method. The resulting mutant is appended to
:math:`P_\mathrm{o}`. In the case of a reproduction, one individual is
selected at random from :math:`P_\mathrm{p}`, cloned and appended to
:math:`P_\mathrm{o}`.
This variation is named *Or* because an offspring will never result from
both operations crossover and mutation. The sum of both probabilities
shall be in :math:`[0, 1]`, the reproduction probability is
1 - *cxpb* - *mutpb*.
"""
assert (cxpb + mutpb) <= 1.0, (
"The sum of the crossover and mutation probabilities must be smaller "
"or equal to 1.0.")
offspring = []
for _ in xrange(lambda_):
op_choice = random.random()
if op_choice < cxpb: # Apply crossover
ind1, ind2 = map(toolbox.clone, random.sample(population, 2))
ind1, ind2 = toolbox.mate(ind1, ind2)
del ind1.fitness.values
offspring.append(ind1)
elif op_choice < cxpb + mutpb: # Apply mutation
ind = toolbox.clone(random.choice(population))
ind, = toolbox.mutate(ind)
del ind.fitness.values
offspring.append(ind)
else: # Apply reproduction
offspring.append(random.choice(population))
return offspring
def eaMuPlusLambda(population, toolbox, mu, lambda_, cxpb, mutpb, ngen,
stats=None, halloffame=None, verbose=__debug__):
"""This is the :math:`(\mu + \lambda)` evolutionary algorithm.
:param population: A list of individuals.
:param toolbox: A :class:`~deap.base.Toolbox` that contains the evolution
operators.
:param mu: The number of individuals to select for the next generation.
:param lambda\_: The number of children to produce at each generation.
:param cxpb: The probability that an offspring is produced by crossover.
:param mutpb: The probability that an offspring is produced by mutation.
:param ngen: The number of generation.
:param stats: A :class:`~deap.tools.Statistics` object that is updated
inplace, optional.
:param halloffame: A :class:`~deap.tools.HallOfFame` object that will
contain the best individuals, optional.
:param verbose: Whether or not to log the statistics.
:returns: The final population
:returns: A class:`~deap.tools.Logbook` with the statistics of the
evolution.
The algorithm takes in a population and evolves it in place using the
:func:`varOr` function. It returns the optimized population and a
:class:`~deap.tools.Logbook` with the statistics of the evolution. The
logbook will contain the generation number, the number of evaluations for
each generation and the statistics if a :class:`~deap.tools.Statistics` is
given as argument. The *cxpb* and *mutpb* arguments are passed to the
:func:`varOr` function. The pseudocode goes as follow ::
evaluate(population)
for g in range(ngen):
offspring = varOr(population, toolbox, lambda_, cxpb, mutpb)
evaluate(offspring)
population = select(population + offspring, mu)
First, the individuals having an invalid fitness are evaluated. Second,
the evolutionary loop begins by producing *lambda_* offspring from the
population, the offspring are generated by the :func:`varOr` function. The
offspring are then evaluated and the next generation population is
selected from both the offspring **and** the population. Finally, when
*ngen* generations are done, the algorithm returns a tuple with the final
population and a :class:`~deap.tools.Logbook` of the evolution.
This function expects :meth:`toolbox.mate`, :meth:`toolbox.mutate`,
:meth:`toolbox.select` and :meth:`toolbox.evaluate` aliases to be
registered in the toolbox. This algorithm uses the :func:`varOr`
variation.
"""
logbook = tools.Logbook()
logbook.header = ['gen', 'nevals'] + (stats.fields if stats else [])
# Evaluate the individuals with an invalid fitness
invalid_ind = [ind for ind in population if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
for ind, fit in zip(invalid_ind, fitnesses):
ind.fitness.values = fit
if halloffame is not None:
halloffame.update(population)
record = stats.compile(population) if stats is not None else {}
logbook.record(gen=0, nevals=len(invalid_ind), **record)
if verbose:
print(logbook.stream)
# Begin the generational process
for gen in range(1, ngen + 1):
# Vary the population
offspring = varOr(population, toolbox, lambda_, cxpb, mutpb)
# Evaluate the individuals with an invalid fitness
invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
for ind, fit in zip(invalid_ind, fitnesses):
ind.fitness.values = fit
# Update the hall of fame with the generated individuals
if halloffame is not None:
halloffame.update(offspring)
# Select the next generation population
population[:] = toolbox.select(population + offspring, mu)
best_ind = halloffame.items[0]
# Update the statistics with the new population
record = stats.compile(population) if stats is not None else {}
logbook.record(gen=gen, nevals=len(invalid_ind), **record)
temp = logbook.stream
if verbose:
print(temp)
print(best_ind)
f = open(time_str + "outputfile.dat", "a")
f.write(str(temp) + '\n')
f.write(str(best_ind) + '\n')
f.close() # 关闭文件
if halloffame.items[0].fitness.values[0] < 20:
print("绝对百分比误差小于20%，完成计算！")
print(best_ind)
return population, logbook
return population, logbook
def varAnd(population, toolbox, cxpb, mutpb):
r"""Part of an evolutionary algorithm applying only the variation part
(crossover **and** mutation). The modified individuals have their
fitness invalidated. The individuals are cloned so returned population is
independent of the input population.
:param population: A list of individuals to vary.
:param toolbox: A :class:`~deap.base.Toolbox` that contains the evolution
operators.
:param cxpb: The probability of mating two individuals.
:param mutpb: The probability of mutating an individual.
:returns: A list of varied individuals that are independent of their
parents.
The variation goes as follow. First, the parental population
:math:`P_\mathrm{p}` is duplicated using the :meth:`toolbox.clone` method
and the result is put into the offspring population :math:`P_\mathrm{o}`. A
first loop over :math:`P_\mathrm{o}` is executed to mate pairs of
consecutive individuals. According to the crossover probability *cxpb*, the
individuals :math:`\mathbf{x}_i` and :math:`\mathbf{x}_{i+1}` are mated
using the :meth:`toolbox.mate` method. The resulting children
:math:`\mathbf{y}_i` and :math:`\mathbf{y}_{i+1}` replace their respective
parents in :math:`P_\mathrm{o}`. A second loop over the resulting
:math:`P_\mathrm{o}` is executed to mutate every individual with a
probability *mutpb*. When an individual is mutated it replaces its not
mutated version in :math:`P_\mathrm{o}`. The resulting :math:`P_\mathrm{o}` is returned.
This variation is named *And* because of its propensity to apply both
crossover and mutation on the individuals. Note that both operators are
not applied systematically, the resulting individuals can be generated from
crossover only, mutation only, crossover and mutation, and reproduction
according to the given probabilities. Both probabilities should be in
:math:`[0, 1]`.
"""
offspring = [toolbox.clone(ind) for ind in population]
# Apply crossover and mutation on the offspring
for i in range(1, len(offspring), 2):
if random.random() < cxpb:
offspring[i - 1], offspring[i] = toolbox.mate(offspring[i - 1],
offspring[i])
del offspring[i - 1].fitness.values, offspring[i].fitness.values
for i in range(len(offspring)):
if random.random() < mutpb:
offspring[i], = toolbox.mutate(offspring[i])
del offspring[i].fitness.values
return offspring
def eaSimple(population, toolbox, cxpb, mutpb, ngen, stats=None,
halloffame=None, verbose=__debug__):
"""This algorithm reproduce the simplest evolutionary algorithm as
presented in chapter 7 of [Back2000]_.
:param population: A list of individuals.
:param toolbox: A :class:`~deap.base.Toolbox` that contains the evolution
operators.
:param cxpb: The probability of mating two individuals.
:param mutpb: The probability of mutating an individual.
:param ngen: The number of generation.
:param stats: A :class:`~deap.tools.Statistics` object that is updated
inplace, optional.
:param halloffame: A :class:`~deap.tools.HallOfFame` object that will
contain the best individuals, optional.
:param verbose: Whether or not to log the statistics.
:returns: The final population
:returns: A class:`~deap.tools.Logbook` with the statistics of the
evolution
The algorithm takes in a population and evolves it in place using the
:meth:`varAnd` method. It returns the optimized population and a
:class:`~deap.tools.Logbook` with the statistics of the evolution. The
logbook will contain the generation number, the number of evaluations for
each generation and the statistics if a :class:`~deap.tools.Statistics` is
given as argument. The *cxpb* and *mutpb* arguments are passed to the
:func:`varAnd` function. The pseudocode goes as follow ::
evaluate(population)
for g in range(ngen):
population = select(population, len(population))
offspring = varAnd(population, toolbox, cxpb, mutpb)
evaluate(offspring)
population = offspring
As stated in the pseudocode above, the algorithm goes as follow. First, it
evaluates the individuals with an invalid fitness. Second, it enters the
generational loop where the selection procedure is applied to entirely
replace the parental population. The 1:1 replacement ratio of this
algorithm **requires** the selection procedure to be stochastic and to
select multiple times the same individual, for example,
:func:`~deap.tools.selTournament` and :func:`~deap.tools.selRoulette`.
Third, it applies the :func:`varAnd` function to produce the next
generation population. Fourth, it evaluates the new individuals and
compute the statistics on this population. Finally, when *ngen*
generations are done, the algorithm returns a tuple with the final
population and a :class:`~deap.tools.Logbook` of the evolution.
.. note::
Using a non-stochastic selection method will result in no selection as
the operator selects *n* individuals from a pool of *n*.
This function expects the :meth:`toolbox.mate`, :meth:`toolbox.mutate`,
:meth:`toolbox.select` and :meth:`toolbox.evaluate` aliases to be
registered in the toolbox.
.. [Back2000] Back, Fogel and Michalewicz, "Evolutionary Computation 1 :
Basic Algorithms and Operators", 2000.
"""
logbook = tools.Logbook()
logbook.header = ['gen', 'nevals'] + (stats.fields if stats else [])
# Evaluate the individuals with an invalid fitness
invalid_ind = [ind for ind in population if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
for ind, fit in zip(invalid_ind, fitnesses):
ind.fitness.values = fit
if halloffame is not None:
halloffame.update(population)
record = stats.compile(population) if stats else {}
logbook.record(gen=0, nevals=len(invalid_ind), **record)
if verbose:
print(logbook.stream)
# Begin the generational process
for gen in range(1, ngen + 1):
# Select the next generation individuals
offspring = toolbox.select(population, len(population))
# Vary the pool of individuals
offspring = varAnd(offspring, toolbox, cxpb, mutpb)
# Evaluate the individuals with an invalid fitness
invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
for ind, fit in zip(invalid_ind, fitnesses):
ind.fitness.values = fit
# Update the hall of fame with the generated individuals
if halloffame is not None:
halloffame.update(offspring)
# Replace the current population by the offspring
population[:] = offspring
# Append the current generation statistics to the logbook
best_ind = halloffame.items[0]
record = stats.compile(population) if stats else {}
logbook.record(gen=gen, nevals=len(invalid_ind), **record)
temp = logbook.stream
if verbose:
print(temp)
print(best_ind)
f = open(time_str + "outputfile.dat", "a")
f.write(str(temp) + '\n')
f.write(str(best_ind) + '\n')
f.close() # 关闭文件
if halloffame.items[0].fitness.values[0] < 20:
print("绝对百分比误差小于20%，完成计算！")
print(best_ind)
return population, logbook
return population, logbook
def main():
random.seed(318)
### 参数设置
cxpb = 0.5 # 交叉概率
mutpb = 0.15 # 变异概率0
ngen = 1500 # 迭代次数
popSize = 150 # 种群规模
### 生成初始种群 - 使用新的验证方法
print("创建初始种群...")
pop = create_valid_population(popSize)
print("初始种群创建完成")
# 验证初始种群中的所有个体
for i, ind in enumerate(pop):
if not check_pow_arguments(ind):
print(f"警告: 初始种群中的个体 {i} 无效!")
print(str(ind))
hof = tools.HallOfFame(1)
invalid_fit = [ind for ind in pop if not ind.fitness.valid]
fitnesses = tool.map(tool.evaluate, invalid_fit) # 求初始种群的每个个体的适应度值，是一个list
for ind, fit in zip(invalid_fit, fitnesses):
ind.fitness.values = fit # 给每个ind的fitness赋值
hof.update(pop)
best_inds = [] # 记录每一代的最优个体
best_ind = hof.items[0]
best_inds.append(best_ind)
stats = tools.Statistics(lambda ind: ind.fitness.values)
stats.register('avg', numpy.mean, axis=0)
### axis = 0表示数组的列求平均；axis = 1 表示数组的行求平均
stats.register('std', numpy.std, axis=0)
stats.register('min', numpy.min, axis=0)
stats.register('max', numpy.max, axis=0)
### record是一个字典：record{‘avg’：[],'std':[],...},每个value是一个list
record = stats.compile(pop)
for key in record:
record[key] = record[key][0]
logbook = tools.Logbook()
logbook.record(gen=0, eval=popSize, best_ind=best_ind, **record)
logbook.header = 'gen', 'min', 'max', 'avg', 'std', 'eval', 'best_ind'
print(logbook)
print('--------开始迭代--------')
for g in range(ngen):
### select 选出popSize个个体进入下一代
# print('The', g + 1, 'step')
offspring = tool.select(pop, len(pop))
# offspring = list(map(tool.clone, offspring))
# offspring = tool.map(tool.clone, offspring)
offspring = [tool.clone(ind) for ind in offspring]
# Apply crossover and mutation on the offspring
for i in range(1, len(offspring), 2):
if random.random() < cxpb:
offspring[i - 1], offspring[i] = tool.mate(offspring[i - 1], offspring[i])
del offspring[i - 1].fitness.values, offspring[i].fitness.values
for i in range(len(offspring)):
if random.random() < mutpb:
offspring[i], = tool.mutate(offspring[i])
del offspring[i].fitness.values
# Evaluate the individuals with an invalid fitness
invalid_fit = [ind for ind in offspring if not ind.fitness.valid]
fitnesses = tool.map(tool.evaluate, invalid_fit)
for ind, fit in zip(invalid_fit, fitnesses):
ind.fitness.values = fit
hof.update(offspring) # 找到本代的最优个体
pop[:] = offspring # 更新种群
best_ind = hof.items[0]
best_inds.append(best_ind)
record = stats.compile(pop) # 数据统计
for key in record:
record[key] = record[key][0]
logbook.record(gen=g + 1, eval=len(invalid_fit), best_ind=best_ind, **record)
temp = logbook.stream
print(temp)
# print logbook.stream
ff = open(time_str + "outputfile.dat", "a")
ff.write(str(temp) + '\n')
# f.write(str(best_ind) + '\n')
ff.close() # 关闭文件
if hof.items[0].fitness.values[0] < 20:
print("绝对百分比误差小于20%，完成计算！")
print(best_ind)
print('--------迭代结束-------')
break
print('--------迭代结束-------')
# print(logbook)
# pop = tool.population(n=150)
# hof = tools.HallOfFame(1)
#
# stats_fit = tools.Statistics(lambda ind: ind.fitness.values)
# stats_size = tools.Statistics(len)
# mstats = tools.MultiStatistics(fitness=stats_fit, size=stats_size)
# mstats.register("avg", numpy.mean)
# mstats.register("std", numpy.std)
# mstats.register("min", numpy.min)
# mstats.register("max", numpy.max)
# pop, log = eaSimple(pop, tool, 0.5, 0.1, 5000, stats=mstats, halloffame=hof, verbose=True)
# pop, log = eaMuPlusLambda(pop, tool, 20, 40, 0.5, 0.1, 5000,
# stats=mstats, halloffame=hof, verbose=True)
# print log
# return pop, log, hof
drawBestInd(expr=hof.items[0])
# 画出实际函数和优化的函数
# temp = numpy.zeros((len(data)))
temp = numpy.zeros((len(data), 5))
func = tool.compile(expr=hof.items[0])
for i in range(len(data)):
temp[i, 1] = (func(Fr[i], We[i], Re[i], Bo[i]) + 1) * dPa[
i] # func(FrG, Bo, Rel, ReG, dPa) func(Fr, Bo, f, Pred, x, N)
x1 = numpy.arange(0, (max(data[:, 65])), 10)
y1 = 1.3 * x1
y2 = 0.7 * x1
# func = tool.compile(expr=hof.items[0])
# z = [func(x) for x in points]
temp[:, 0] = data[:, 65]
temp[:, 2] = data[:, 66]
# temp[:, 3] = data[:, 7]
# temp[:, 4] = data[:, 8]
filename = time_str + '_data.dat'
text_save(filename, temp)
# f1 = open(time_str + "outputfile.dat", "a")
# f1.write( + '\n')
# f1.close
ax1 = plt.subplot(211)
plt.plot(data[:, 65], temp[:, 1], 'o', markeredgecolor='b', markerfacecolor='none', label='EXP-PRE-ENV')
plt.plot(x1, y1, '--', color='k') # , label='+30%')
plt.plot(x1, y2, '--', color='k') # , label='-30%')
plt.legend(loc='best')
plt.xscale('log')
plt.yscale('log')
ax1 = plt.subplot(212)
plt.plot(data[:, 65], data[:, 66], 'o', markeredgecolor='r', markerfacecolor='none', label='EXP-PRE_ANA')
plt.plot(x1, y1, '--', color='k') # , label='+30%')
plt.plot(x1, y2, '--', color='k') # , label='-30%')
plt.legend(loc='best')
plt.xscale('log')
plt.yscale('log')
# ax1 = plt.subplot(223)
# plt.plot(data[:, 5], data[:, 7], 'o', markeredgecolor='c', markerfacecolor='none', label='EXP-PRE_MSH')
# plt.plot(x1, y1, '--', color='k')#, label='+30%')
# plt.plot(x1, y2, '--', color='k')#, label='-30%')
# plt.legend(loc='best')
# plt.xscale('log')
# plt.yscale('log')
# ax1 = plt.subplot(224)
# plt.plot(data[:, 5], data[:, 8], 'o', markeredgecolor='g', markerfacecolor='none', label='EXP-PRE_M&H')
plt.plot(x1, y1, '--', color='k') # , label='+30%')
plt.plot(x1, y2, '--', color='k') # , label='-30%')
# plt.legend(loc='best')
# plt.xscale('log')
# plt.yscale('log')
# ax.set_aspect(1)
# curr_time = datetime.datetime.now()
# time_str=datetime.datetime.strftime(curr_time, '%Y%m%d%H%M')
plt.savefig(time_str + 'res.pdf')
plt.show()
main()
最新发布
介绍了一款能够在大量文件中轻松搜索并替换文本的工具，对于需要在众多文件中进行文本修改的任务非常实用。
扫一扫
743
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
