From b95a55d43d172c0b07f621aa2836d9be56b38a54 Mon Sep 17 00:00:00 2001 From: XinWang2021 <9833391+xinwang2021@user.noreply.gitee.com> Date: Wed, 1 Dec 2021 12:42:27 +0000 Subject: [PATCH 1/2] =?UTF-8?q?=E6=96=B0=E5=BB=BA=20=E7=8E=8B=E9=91=AB-?= =?UTF-8?q?=E5=9F=BA=E4=BA=8E=E5=A4=9A=E5=B1=82=E7=A5=9E=E7=BB=8F=E7=BD=91?= =?UTF-8?q?=E7=BB=9C=E7=9A=84=E5=90=B8=E6=94=B6=E8=BE=B9=E7=95=8C=E6=9D=A1?= =?UTF-8?q?=E4=BB=B6=E4=B8=8B=E7=9A=84=E4=B8=80=E7=BB=B4=E8=96=9B=E5=AE=9A?= =?UTF-8?q?=E8=B0=94=E6=96=B9=E7=A8=8B=E6=95=B0=E5=80=BC=E8=A7=A3=E6=B3=95?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- .../.keep" | 0 1 file changed, 0 insertions(+), 0 deletions(-) create mode 100644 "code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/.keep" diff --git "a/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/.keep" "b/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/.keep" new file mode 100644 index 0000000..e69de29 -- Gitee From 3ae9869001263eee29922cfd0283ec1066218ca4 Mon Sep 17 00:00:00 2001 From: XinWang2021 <9833391+xinwang2021@user.noreply.gitee.com> Date: Wed, 1 Dec 2021 12:43:14 +0000 Subject: [PATCH 2/2] =?UTF-8?q?=E8=96=9B=E5=AE=9A=E8=B0=94=E6=96=B9?= =?UTF-8?q?=E7=A8=8B=E6=95=B0=E5=80=BC=E6=B1=82=E8=A7=A3?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- .../Multi_layers_net.py" | 39 ++ .../dataset.py" | 107 +++++ .../train_model.py" | 384 ++++++++++++++++++ 3 files changed, 530 insertions(+) create mode 100644 "code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/Multi_layers_net.py" create mode 100644 "code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/dataset.py" create mode 100644 "code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/train_model.py" diff --git "a/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/Multi_layers_net.py" "b/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/Multi_layers_net.py" new file mode 100644 index 0000000..8d4cf26 --- /dev/null +++ "b/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/Multi_layers_net.py" @@ -0,0 +1,39 @@ +# author: Xin Wang 王鑫 +# email: 2021202010045@whu.edu.cn +# 2021-11 +# coding: utf-8 + +import tensorflow as tf +import numpy as np + + +def add_layer2(x, input_dim=1, output_dim=1, name_scope='hidden'): + with tf.compat.v1.variable_scope(name_scope, reuse=tf.compat.v1.AUTO_REUSE): + ua_w = tf.compat.v1.get_variable(name='ua_w', shape=(input_dim, output_dim), dtype=tf.float32, + initializer=tf.compat.v1.glorot_uniform_initializer()) # Xavier初始化相同 + ua_b = tf.compat.v1.get_variable(name='ua_b', shape=(output_dim,), dtype=tf.float32, + initializer=tf.compat.v1.glorot_uniform_initializer()) + z = tf.matmul(x, ua_w) + ua_b + output_z = tf.nn.tanh(z) + + return output_z + + +def univAprox2(x0, hidden_units=[10, 20, 40], input_dim=1, output_dim_final=1): + hidden_num = len(hidden_units) + add_hidden = [input_dim] + hidden_units + output = x0 + for i in range(hidden_num): + input_dim = add_hidden[i] + output_dim = add_hidden[i + 1] + print('input_dim:%s, output_dim:%s' % (input_dim, output_dim)) + name_scope = 'hidden' + np.str(i + 1) + output = add_layer2(output, input_dim, output_dim, name_scope=name_scope) + + ua_we = tf.compat.v1.get_variable(name='ua_we', shape=(add_hidden[-1], output_dim_final), dtype=tf.float32, + initializer=tf.compat.v1.glorot_uniform_initializer()) + ua_be = tf.compat.v1.get_variable(name='ua_be', shape=(output_dim_final,), dtype=tf.float32, + initializer=tf.compat.v1.glorot_uniform_initializer()) + + z = tf.matmul(output, ua_we) + ua_be + return z diff --git "a/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/dataset.py" "b/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/dataset.py" new file mode 100644 index 0000000..baa9ff9 --- /dev/null +++ "b/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/dataset.py" @@ -0,0 +1,107 @@ +# author: Xin Wang 王鑫 +# email: 2021202010045@whu.edu.cn +# 2021-11 +# coding: utf-8 + +import numpy as np + +Sig = 0.04 +K = 2 + + +# 随机取值，返回自变量矩阵[t, x], size: (n_size, 1 + x_dim) +def get_rand_t_and_x(n_size, x_dim, xrange=[0, 1], trange=[0, 1]): + xxs = np.random.rand(n_size, x_dim) * (xrange[1] - xrange[0]) + xrange[0] + tts = np.random.rand(n_size, 1) * (trange[1] - trange[0]) + trange[0] + + return np.c_[tts, xxs] + + +# t = 0时，边界点随机剖分取值，返回t=0是的边界点自变量数据矩阵[t=0，x], size: (n_size, 1 + x_dim) +def get_t_boundary(n_size, x_dim, xrange=[0, 1], trange=[0, 1]): + t_xx_b = get_rand_t_and_x(n_size, x_dim, xrange, trange) + t_xx_b[:, 0] = trange[0] + + return t_xx_b + + +# 空间上边界点随机剖分取值，返回空间边界点自变量数据矩阵[t, x = xrange], size: (n_size * x_dim * 2, 1 + x_dim) +def get_x_boundary(n_size, x_dim, xrange=[0, 1], trange=[0, 1]): + t_xx_b = get_rand_t_and_x(n_size * x_dim * 2, x_dim, xrange, trange) + for ii in range(x_dim): + t_xx_b[2 * ii * n_size:(2 * ii + 1) * n_size, ii + 1] = xrange[0] + t_xx_b[(2 * ii + 1) * n_size:2 * (ii + 1) * n_size, ii + 1] = xrange[1] + + return t_xx_b + + +# 设置目标函数，利用已知的函数（利用其他偏微分方程解法得到）求解得到的监督数据y = [y.real, y.imag]，size: (length(xx_in), output_dim) +# 注意：由于没有合适的方式处理复数域上的神经网络，选择将复数的实部和虚部分开计算，不影响具体的计算结果。 +def get_target_func(xx_in): + t_tmp = xx_in[:, 0].reshape(-1, 1) + x_tmp = xx_in[:, 1:] + if x_tmp.shape[1] == 1: + k = K + sigma = Sig + tmp = 1j * k * (x_tmp - k * t_tmp) - (x_tmp - 2 * k * t_tmp) ** 2 / (4 * (sigma + 1j * t_tmp)) + y_tmp = np.exp( + 1j * k * (x_tmp - k * t_tmp) - (x_tmp - 2 * k * t_tmp) ** 2 / (4 * (sigma + 1j * t_tmp))) / np.sqrt( + sigma + 1j * t_tmp) + else: + y_tmp = 0 + y_tmp += np.sum(xx_in ** 2, axis=1, keepdims=True) + + return np.c_[y_tmp.real, y_tmp.imag] + + +# 设置空间项边界函数，根据边界条件求解(t, x = xrange)时的监督数据y = [y.real, y.imag]，size: (length(xx_boundary), output_dim) +def get_y_boundary(xx_boundary): + y_tmp = np.zeros([xx_boundary.shape[0], 1]) + + return np.c_[y_tmp.real, y_tmp.imag] + + +# 设置关于空间项的二阶偏微分计算函数，利用已知的函数（利用其他偏微分方程解法得到）求解得到的监督数据grad2y = [y.real, y.imag]，size: (length(xx_in), output_dim) +def get_grad2_x_target(xx_in): + t_tmp = xx_in[:, 0].reshape(-1, 1) + x_tmp = xx_in[:, 1:] + k = K + sigma = Sig + u_tmp = get_target_func(xx_in)[:, 0].reshape(-1, 1) + 1j * get_target_func(xx_in)[:, 1].reshape(-1, 1) + # print(u_tmp) + grad2y_x_tmp = (1j * k - (x_tmp - 2 * k * t_tmp) / (2 * (sigma + 1j * t_tmp))) ** 2 * u_tmp - u_tmp / ( + 2 * (sigma + 1j * t_tmp)) + + return np.c_[grad2y_x_tmp.real, grad2y_x_tmp.imag] + + +# 设置关于时间项的微分计算函数，利用已知的函数（利用其他偏微分方程解法得到）求解得到的监督数据grady = [y.real, y.imag]，size: (length(xx_in), output_dim) +def get_grad_t_target(xx_in): + t_tmp = xx_in[:, 0].reshape(-1, 1) + x_tmp = xx_in[:, 1:] + # print(t_tmp.shape, x_tmp) + k = K + sigma = Sig + u_tmp = get_target_func(xx_in)[:, 0].reshape(-1, 1) + 1j * get_target_func(xx_in)[:, 1].reshape(-1, 1) + grady_t_tmp = (k * (x_tmp - 2 * k * t_tmp) / (sigma + 1j * t_tmp) + 1j * (x_tmp - 2 * k * t_tmp) ** 2 / ( + 4 * (sigma + 1j * t_tmp) ** 2) - 1j * k ** 2) * u_tmp - 1j * u_tmp / ( + 2 * (sigma + 1j * t_tmp)) + + return np.c_[grady_t_tmp.real, grady_t_tmp.imag] + + +# 获取实验数据集 +def get_xy(n_size_in, n_size_b_x, n_size_b_t, x_dim, x_range, t_range): + t_and_x_inside = get_rand_t_and_x(n_size_in, x_dim, xrange=x_range, trange=t_range) + y_inside = get_target_func(t_and_x_inside) + + xx_b = get_x_boundary(n_size_b_x, x_dim, xrange=x_range, trange=t_range) + yy_b_x = get_y_boundary(xx_b) + + tt_b = get_t_boundary(n_size_b_t, x_dim, xrange=x_range, trange=t_range) + yy_b_t = get_target_func(tt_b) + + g2d_y_x_in = get_grad2_x_target(t_and_x_inside) + gd_y_t_in = get_grad_t_target(t_and_x_inside) + + return t_and_x_inside, y_inside, xx_b, yy_b_x, tt_b, yy_b_t, g2d_y_x_in, gd_y_t_in diff --git "a/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/train_model.py" "b/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/train_model.py" new file mode 100644 index 0000000..2022c3a --- /dev/null +++ "b/code/2021_autumn/\347\216\213\351\221\253-\345\237\272\344\272\216\345\244\232\345\261\202\347\245\236\347\273\217\347\275\221\347\273\234\347\232\204\345\220\270\346\224\266\350\276\271\347\225\214\346\235\241\344\273\266\344\270\213\347\232\204\344\270\200\347\273\264\350\226\233\345\256\232\350\260\224\346\226\271\347\250\213\346\225\260\345\200\274\350\247\243\346\263\225/train_model.py" @@ -0,0 +1,384 @@ +# author: Xin Wang 王鑫 +# email: 2021202010045@whu.edu.cn +# 2021-11-7 +# coding: utf-8 + + +import os + +import pickle +import time +import shutil +import tensorflow as tf +import numpy as np +import matplotlib.pyplot as plt +from Sch_by_Wang.dataset import * +from Sch_by_Wang.Multi_layers_net import * + +from mpl_toolkits.mplot3d import Axes3D +from matplotlib import cm +from matplotlib.ticker import LinearLocator, FormatStrFormatter + +Leftp = 0.18 +Bottomp = 0.18 +Widthp = 0.88 - Leftp +Heightp = 0.9 - Bottomp +pos = [Leftp, Bottomp, Widthp, Heightp] + + +def mkdir(fn): + if not os.path.isdir(fn): + os.mkdir(fn) + + +def mySaveFig(pltm, fntmp, fp=0, ax=0, isax=0, iseps=0, isShowPic=0): + if isax == 1: + pltm.rc('xtick', labelsize=18) + pltm.rc('ytick', labelsize=18) + ax.set_position(pos, which='both') + fnm = '%s.png' % (fntmp) + pltm.savefig(fnm) + if iseps: + fnm = '%s.eps' % (fntmp) + pltm.savefig(fnm, format='eps', dpi=600) + if fp != 0: + fp.savefig("%s.pdf" % (fntmp), bbox_inches='tight') + if isShowPic == 1: + pltm.show() + elif isShowPic == -1: + return + else: + pltm.close() + + +R = {} ### used for saved all parameters and data + +### mkdir a folder to save all output +BaseDir = 'result/' + +subFolderName = '%s' % (int(np.absolute(np.random.normal([1]) * 100000)) // int(1)) +# np.absolute() Calculate the absolute value element-wise.逐个计算每个元素的模值 +FolderName = '%s%s/' % (BaseDir, subFolderName) +mkdir(BaseDir) +mkdir(FolderName) +mkdir('%smodel/' % (FolderName)) +print(FolderName) + +R['FolderName'] = FolderName +R['input_dim'] = 2 +R['x_dim'] = 1 +R['output_dim'] = 2 +R['hidden_units'] = [200, 100, 20] +R['learning_rate'] = 2e-5 +R['learning_rateDecay'] = 5e-7 # 学习率衰减 +R['issavemodel'] = True + +plotepoch = 500 +R['X_in_train_size'] = 2000 ### training size of inside value +R['x_b_train_size'] = 500 ### training size of x boundary +R['t_b_train_size'] = 500 ### training size of t boundary + +R['X_in_test_size'] = 500 ### test size of inside value +R['x_b_test_size'] = 100 ### test size of x boundary +R['t_b_test_size'] = 100 ### test size of t boundary + +R['beta'] = 10 +R['gamma'] = 10 + +# 真解参数设置 +R['sigma'] = 0.04 +R['k'] = 2 + +# 自变量取值范围 +R['xrange'] = [-3, 3] +R['trange'] = [0, 2] +R['epsilon'] = 1 + +R['tol'] = 1e-4 +R['Total_Step'] = 600000 ### the training step. Set a big number, if it converges, can manually stop training + +print('the test parameter is ', R) + + +def get_sample_train(): + R['train_t_and_x_inside'], R['train_y_inside'], R['train_xx_b'], R['train_yy_b_x'], R['train_tt_b'], R[ + 'train_yy_b_t'], R['train_g2d_yy_x_in'], R['train_gd_yy_t_in'] = get_xy(R['X_in_train_size'], + R['x_b_train_size'], + R['t_b_train_size'], R['x_dim'], + R['xrange'], R['trange']) + + +def get_sample_test(): + R['test_t_and_x_inside'], R['test_y_inside'], R['test_xx_b'], R['test_yy_b_x'], R['test_tt_b'], R['test_yy_b_t'], R[ + 'test_g2d_yy_x_in'], R['test_gd_yy_t_in'] = get_xy(R['X_in_test_size'], R['x_b_test_size'], R['t_b_test_size'], + R['x_dim'], R['xrange'], R['trange']) + + +def get_sample(): + R['train_t_and_x_inside'], R['train_y_inside'], R['train_xx_b'], R['train_yy_b_x'], R['train_tt_b'], R[ + 'train_yy_b_t'], R['train_g2d_yy_x_in'], R['train_gd_yy_t_in'] = get_xy(R['x_in_train_size'], + R['x_b_train_size'], + R['t_b_train_size']) + Y = np.r_[R['train_y_inside'], R['train_yy_b_x'], R['train_yy_b_t']] + R['Y_min'] = np.min(np.sum(np.square(Y), axis=1)) + R['Y_max'] = np.max(np.sum(np.square(Y), axis=1)) + # print(np.where(R['train_y_inside']==0)) + # print('the shape of all data is :') + # print(R['train_t_and_x_inside'].dtype, R['train_y_inside'].dtype, R['train_xx_b'].shape, R['train_yy_b_x'].shape, + # R['train_tt_b'].shape, R['train_yy_b_t'].shape) + R['test_x_inside'], R['test_y_inside'], R['test_xx_b'], R['test_yy_b_x'], R['test_tt_b'], R['test_yy_b_t'], R[ + 'test_g2d_yy_x_in'], R['test_gd_yy_t_in'] = get_xy(R['x_in_test_size'], R['x_b_test_size'], R['t_b_test_size']) + + +with tf.compat.v1.variable_scope('Graph', reuse=tf.compat.v1.AUTO_REUSE) as scope: + # Our inputs will be a batch of values taken by our functions + # x is the interior points + + tf.compat.v1.disable_eager_execution() + x = tf.compat.v1.placeholder(tf.float32, shape=[None, R['input_dim']], name="x") + # x_len is the length of x + x_len = tf.compat.v1.placeholder_with_default(input=1, shape=[], name='xlen') + # xb is the sample on boundary + xb = tf.compat.v1.placeholder(tf.float32, shape=[None, R['input_dim']], name="xb") + # tb is the sample on boundary + tb = tf.compat.v1.placeholder(tf.float32, shape=[None, R['input_dim']], name="tb") + # b_len is the length of x + b_len = tf.compat.v1.placeholder_with_default(input=1, shape=[], name='blen') + # y_true is the label of interior samples + y_true = tf.compat.v1.placeholder(tf.float32, shape=[None, R['output_dim']], name="y") + # y_true_b is the label of the sample on boundary + y_true_b_x = tf.compat.v1.placeholder(tf.float32, shape=[None, R['output_dim']], name="yb_x") + # y_true_b is the label of the sample on boundary + y_true_b_t = tf.compat.v1.placeholder(tf.float32, shape=[None, R['output_dim']], name="yb_t") + # tg2u_x is the true 2nd x gradient of interior samples + tg2u_x = tf.compat.v1.placeholder(tf.float32, shape=[None, R['output_dim']], name="tg2u_x") + # tgu_t is the true 1nd t gradient of interior samples + tgu_t = tf.compat.v1.placeholder(tf.float32, shape=[None, R['output_dim']], name="tgu_t") + # beta is the penality weight of boundary + beta = tf.compat.v1.placeholder_with_default(input=1.0, shape=[], name='beta') + # _lr is the learning rate of each step. + _lr = tf.compat.v1.placeholder_with_default(input=1e-3, shape=[], name='lr') + # Laplacian + # To compute laplacian, we need to unstack of x, such that + # each dimension is one variable + _x = tf.unstack(x, axis=1) + # unstack each dimension, so each dim is a variable + # We pack all dim together so we can get their predict values + # first, we need to expand each dim to shape of [-1,1] + # then, pack them as the final shape as [-1,1] + _x_ = [tf.expand_dims(tmp, axis=1) for tmp in _x] + _x2 = tf.transpose(tf.stack(_x_))[0] + # all inputs combines both the interior and boundary samples. + _x_0 = tf.concat([_x2, xb, tb], axis=0) + _y = tf.concat([y_true, y_true_b_x, y_true_b_t], axis=0) + # univAprox2 is the DNN that returns predicted labels + y = univAprox2(_x_0, R['hidden_units'], + input_dim=R['input_dim'], + output_dim_final=R['output_dim']) + + # get the 2nd derivative of x + # the real part + y_real = y[:, 0:1] + _grady_x_real = [tf.squeeze(tf.gradients(y_real[0:x_len], tmp)) for tmp in _x] + grad2y_x_real = 0 + for ii in range(R['x_dim']): + # take out each dimension and do grad + grad2y_x_real += tf.stack(tf.gradients(_grady_x_real[ii + 1], _x_[ii + 1]))[0] + # the imaginary part + y_imag = y[:, 1:] + _grady_x_imag = [tf.squeeze(tf.gradients(y_imag[0:x_len], tmp)) for tmp in _x] + grad2y_x_imag = 0 + for ii in range(R['x_dim']): + # take out each dimension and do grad + grad2y_x_imag += tf.stack(tf.gradients(_grady_x_imag[ii + 1], _x_[ii + 1]))[0] + + grad2y = tf.concat([grad2y_x_real, grad2y_x_imag], axis=1) + + # get the first derivative of t . + y_real = y[:, 0:1] + _grady_t_real = tf.transpose(tf.gradients(y_real[0:x_len], _x[0])) + + y_imag = y[:, 1:] + _grady_t_imag = tf.transpose(tf.gradients(y_imag[0:x_len], _x[0])) + _grady_t = tf.concat([_grady_t_real, _grady_t_imag], axis=1) + + loss_t = tf.reduce_mean(tf.reduce_sum(tf.square(_grady_t - tgu_t), axis=1)) + loss_x = tf.reduce_mean(tf.reduce_sum(tf.square(grad2y - tg2u_x), axis=1)) + + lossb = tf.reduce_mean(tf.reduce_sum(tf.square(_y[b_len:] - y[b_len:]), axis=1)) + loss = loss_t + loss_x + beta * lossb + lossy_in = tf.reduce_mean(tf.reduce_sum(tf.square(_y[0:x_len] - y[0:x_len]), axis=1)) + lossy = tf.reduce_mean(tf.reduce_sum(tf.square(_y - y), axis=1)) + + # define train operation using the Adam optimizer + adam = tf.compat.v1.train.AdamOptimizer(learning_rate=_lr) + train_op = adam.minimize(loss) + +config = tf.compat.v1.ConfigProto() +config.gpu_options.allow_growth = True +sess = tf.compat.v1.Session(config=config) + +sess.run(tf.compat.v1.global_variables_initializer()) +saver = tf.compat.v1.train.Saver() + +t0 = time.time() + + +class model(): + def __init__(self): + R['y_train'] = [] + R['y_test'] = [] + R['loss_train'] = [] + R['lossy_train'] = [] + R['lossy_in_rela_train'] = [] + R['losst_train'] = [] + R['lossx_train'] = [] + R['lossb_train'] = [] + R['lossy_in_train'] = [] + if R['issavemodel']: + nametmp = '%smodel/' % (FolderName) + mkdir(nametmp) + saver.save(sess, "%smodel.ckpt" % (nametmp)) + + def run_onestep(self): + get_sample_train() + R['y_train'], losst, lossx, loss_b, lossy_in_tmp, lossy_train_tmp, loss_train_tmp, g1, g2 = sess.run( + [y, loss_t, loss_x, lossb, lossy_in, lossy, loss, _grady_t, grad2y], + feed_dict={x: R['train_t_and_x_inside'], x_len: np.shape(R['train_t_and_x_inside'])[0], xb: R['train_xx_b'], + y_true: R['train_y_inside'], y_true_b_x: R['train_yy_b_x'], tb: R['train_tt_b'], + b_len: np.shape(R['train_tt_b'])[0], + y_true_b_t: R['train_yy_b_t'], tg2u_x: R['train_g2d_yy_x_in'], tgu_t: R['train_gd_yy_t_in'], + beta: R['beta']}) + R['loss_train'].append(loss_train_tmp) + R['lossy_train'].append(lossy_train_tmp) + R['losst_train'].append(losst) + R['lossx_train'].append(lossx) + R['lossb_train'].append(loss_b) + R['lossy_in_train'].append(lossy_in_tmp) + _ = sess.run(train_op, feed_dict={x: R['train_t_and_x_inside'], x_len: np.shape(R['train_t_and_x_inside'])[0], + xb: R['train_xx_b'], y_true: R['train_y_inside'], + y_true_b_x: R['train_yy_b_x'], tb: R['train_tt_b'], + b_len: np.shape(R['train_tt_b'])[0], + y_true_b_t: R['train_yy_b_t'], tg2u_x: R['train_g2d_yy_x_in'], + tgu_t: R['train_gd_yy_t_in'], beta: R['beta'], _lr: R['learning_rate']}) + R['learning_rate'] = R['learning_rate'] * (1 - R['learning_rateDecay']) + + def run(self, step_n=1): + if R['issavemodel']: + nametmp = '%smodel/model.ckpt' % (FolderName) + saver.restore(sess, nametmp) + for ii in range(step_n): + self.run_onestep() + if R['loss_train'][-1] < R['tol']: + print('model end, error:%s' % (R['lossy_train'][-1])) + self.plotloss() + self.savefile() + + if R['issavemodel']: + nametmp = '%smodel/' % (FolderName) + shutil.rmtree(nametmp) + saver.save(sess, "%smodel.ckpt" % (nametmp)) + break + + if ii % plotepoch == 0: + get_sample_test() + R[ + 'y_test'], losst_tmp, lossx_tmp, loss_b_tmp, lossy_in_tmp, lossy_test_tmp, loss_test_tmp, g1, g2 = sess.run( + [y, loss_t, loss_x, lossb, lossy_in, lossy, loss, _grady_t, grad2y], + feed_dict={x: R['test_t_and_x_inside'], x_len: np.shape(R['test_t_and_x_inside'])[0], + xb: R['test_xx_b'], y_true: R['test_y_inside'], y_true_b_x: R['test_yy_b_x'], + tb: R['test_tt_b'], b_len: np.shape(R['test_tt_b'])[0], + y_true_b_t: R['test_yy_b_t'], tg2u_x: R['test_g2d_yy_x_in'], tgu_t: R['test_gd_yy_t_in'], + beta: R['beta']}) + print('time elapse: %.3f' % (time.time() - t0)) + print('epoch: %d, train loss: %.9f | test loss: %.9f' % (ii, R['loss_train'][-1], loss_test_tmp)) + print('epoch: %d, train lossy: %.9f, lossy_in: %.9f | test lossy: %.9f, lossy_in: %.9f' % ( + ii, R['lossy_train'][-1], R['lossy_in_train'][-1], lossy_test_tmp, lossy_in_tmp)) + # print('model, epoch: %d, train losst: %.9f, lossx: %.9f, lossb: %.9f' % (ii, R['losst_train'][-1], R['lossx_train'][-1], R['lossb_train'][-1])) + # print('the predict y range: [(%.9f), (%.9f)], the true y range: [(%.9f), (%.9f)]' % (np.min(np.sum(np.square(R['y_train']), axis=1)), np.max(np.sum(np.square(R['y_train']), axis=1)), R['Y_min'], R['Y_max'])) + if ii == step_n - 1: + tmp_X = np.r_[R['test_t_and_x_inside'], R['test_xx_b'], R['test_tt_b']] + np.savetxt("{}pre.txt".format(R['FolderName']), np.c_[tmp_X, R['y_test']], fmt='%f', delimiter=' ') + np.savetxt("{}y_true.txt".format(R['FolderName']), + np.c_[tmp_X, np.r_[R['test_y_inside'], R['test_yy_b_x'], R['test_yy_b_t']]], fmt='%f', + delimiter=' ') + + self.plotloss() + self.savefile() + if R['issavemodel']: + nametmp = '%smodel/' % (FolderName) + shutil.rmtree(nametmp) + saver.save(sess, "%smodel.ckpt" % (nametmp)) + + def plotloss(self): + plt.figure() + ax = plt.gca() + y2 = R['loss_train'] + plt.plot(y2, 'g*', label='Train') + ax.set_xscale('log') + ax.set_yscale('log') + plt.legend(fontsize=18) + plt.title('loss', fontsize=15) + fntmp = '%sloss' % (FolderName) + mySaveFig(plt, fntmp, ax=ax, isax=1, iseps=0) + + plt.figure() + ax = plt.gca() + y2 = R['lossy_train'] + plt.plot(y2, 'g*', label='Train') + ax.set_xscale('log') + ax.set_yscale('log') + plt.legend(fontsize=18) + plt.title('lossy', fontsize=15) + fntmp = '%slossy' % (FolderName) + mySaveFig(plt, fntmp, ax=ax, isax=1, iseps=0) + + for jj in range(R['output_dim']): + # x: R['test_t_and_x_inside'], x_len: np.shape(R['test_t_and_x_inside'])[0], + # xb: R['test_xx_b'], y_true: R['test_y_inside'], y_true_b_x: R['test_yy_b_x'], + # tb: R['test_tt_b'], b_len: np.shape(R['test_tt_b'])[0], + # y_true_b_t: R['test_yy_b_t'], tg2u_x: R['test_g2d_yy_x_in'], tgu_t: R['test_gd_yy_t_in'], + # beta: R['beta']} + if R['input_dim'] == 1: + plt.figure() + ax = plt.gca() + plt.plot(R['test_t_and_x_inside'], R['y_test'][0:np.shape(R['test_t_and_x_inside'])[0], jj], 'ro', + label='train') + plt.plot(R['test_t_and_x_inside'], R['test_y_inside'][:, jj], 'g*', label='True') + # ax.set_xscale('log') + # ax.set_yscale('log') + plt.legend(fontsize=18) + plt.title('ytest_{}'.format(jj), fontsize=15) + fntmp = '%sy' % (FolderName) + mySaveFig(plt, fntmp, ax=ax, isax=1, iseps=0) + + if R['input_dim'] == 2: + fp = plt.figure() + ax = fp.gca(projection='3d') + ax.scatter(R['test_t_and_x_inside'][:, 0], R['test_t_and_x_inside'][:, 1], R['test_y_inside'][:, jj], + label='true') + ax.scatter(R['test_t_and_x_inside'][:, 0], R['test_t_and_x_inside'][:, 1], + R['y_test'][0:np.shape(R['test_t_and_x_inside'])[0], jj], label='dnn') + # Customize the z axis. + # ax.set_zlim(-2.01, 2.01) + ax.zaxis.set_major_locator(LinearLocator(5)) + ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f')) + plt.legend(fontsize=18) + fntmp = '%sy_%d' % (FolderName, jj) + mySaveFig(plt, fntmp, ax=ax, isax=1, iseps=0) + + def savefile(self): + with open('%s/objs.pkl' % (FolderName), 'wb') as f: # Python 3: open(..., 'wb') + pickle.dump(R, f, protocol=4) + + text_file = open("%s/Output.txt" % (FolderName), "w") + for para in R: + if np.size(R[para]) > 20: + continue + text_file.write('%s: %s\n' % (para, R[para])) + + text_file.close() + + +model1 = model() +model1.run(100001) -- Gitee